Computer Science Experimentation

Saturday, December 25, 2010

WCF with F# interactive – Data Server example

Window Communication Foundation provides a runtime environment for services with the following characteristics:
- Exposes CLR types as services and consumes other services as CLR types. The conversion for the different protocol is done behind the scene, so that all the coding is done as you are in .NET environment.
- Implements the concept of remote objects with concurrency management (introduced by .NET Remoting).
- Implements a set of industrial standard protocols (TCP, HTCP, WSDL, REST)
- Implements fast communication protocol as MSWindows pipes
- Allows multiple bindings for the same .NET service code.
- Provides Hosting.
- Manages Fault, Security and Transactions.

The following example implements a client/server data-server.
Data-server is a service that can be used alone or added to an application engine. It provides efficient storage and retrieve of data. It allows usage of groups of points to improve the reading efficiency. The following example serves float type data. Other types of data as integer, double and string can be added. Another option is to create a generic contract and use multiple servers for each type of data.

The data server provides the same features that OPC-DA, with the addition of the possibility of extra bindings to pipes, html and other web-services protocols. The same way, OPC-DA .NET 3.0 (1.2) is a wrapper on the OPC server (based on COM) to provide the same extra functionalities.

In order to allow an easy environment for development and visualization and dynamic modifications of the application, this example was developed as scripts that are executed in the F# interactive environment. To add intellisense, it is used 2 instances of Visual Studio with F# Interactive.
The first VS is used to develop the service contract and to develop and execute the server in F# interactive.

The second VS is used to develop and execute the client in F# interactive.

This example, a MSWindows pipe binding without any security checks was used.

In order to allow all code development inside F# scripts, the server and client binding configuration was done programmatically. For the same reason, the client works directly with the channel.

The VS were opened in administrator mode to avoid any security blocking (for the file C:\Program Files (x86)\Microsoft Visual Studio 10.0\Common7\IDE\devenv.exe).

The file DS_Contract.fsx implements the data server contract.

This example code includes:
- Definition of DataItem type.
- Definition of IDataServer type with abstract types that includes the signature of all remote methods (interface).
- Definition of DataServer type with the implementation of remote and local methods. The local methods can be only used inside the server host.

This server can be used as stand-alone or embedded in another process, since it executes in a separeted thread to respond to the clients. For the stand-alone case, all the management of data items will be done remotely. For the embedded case, the local application can manage the data items.

The data items are kept in a immutable list, so it can’t be removed or updated after it is inserted. In order to modify it, the application needs to generate the full list of data items. This property provides maximum and consistent data reading speed.

The data server uses a dictionary to relate the item name (tag) to the position of the data list.

In order increase reading speed, the client should create groups of data. The group will hold an array to the data list positions. The group should be deleted if not needed anymore.

DS_Contract.fsx code:
#r "System.ServiceModel"
#r "System.Runtime.Serialization"

//namespace RTE.DataServer

open System
open System.ServiceModel
open System.Collections.Concurrent
open System.Collections.Generic
open System.Runtime.Serialization

    [<DataContract(Namespace="http://RTE/DataServer")>] 
    type DataItem()=
        let mutable Tag_=""
        let mutable Description_=""
        let mutable Value_=0.0
        let mutable Status_=0
        [<datamember>]
        member this.Tag with get()=Tag_ and set x = Tag_<-x
        [<datamember>]
        member this.Description with get()=Description_ and set x = Description_<-x
        [<datamember>]
        member this.Value with get()=Value_ and set x = Value_<-x
        [<datamember>]
        member this.Status with get()=Status_ and set x = Status_<-x

        override this.ToString()=
            sprintf "Tag=%s, Description=%s Value=%f, Status=%i" Tag_ Description_ Value_ Status_ 


    //service contract
    [<Servicecontract (namespace="http://RTE/DataServer" )>]
    type IDataServer=
        [<operationcontract>]
        abstract SetDataValueByName : tag:string -> value:float -> int
        [<operationcontract>]
        abstract SetDataStatusByName : tag:string -> value:int -> int
        [<operationcontract>]
        abstract GetDataValueByName : tag:string -> float*int

        [<operationcontract>]
        abstract RegisterData : tag:DataItem -> int
        [<operationcontract>]
        abstract RegisterDataArray : tags:DataItem[] -> int

        [<operationcontract>]
        abstract AddGroup : group:string -> tags:string[] -> int*int[]
        [<operationcontract>]
        abstract DeleteGroup : group:string -> int
        [<operationcontract>]
        abstract GetGroup : group:string -> int*DataItem[]
        [<operationcontract>]
        abstract GetGroupArray : group:string -> int*float[]*int[]


    //service implementation
    [<Servicebehavior(instancecontextmode=instancecontextmode.single,concurrencymode=concurrencymode.single)>]
    type DataServer()=
        //Constructor & Fields
        let DataDict= new ConcurrentDictionary()
        let DataList= new List()
        let GroupDict= new ConcurrentDictionary()
        let lockObj= ref 0
        //--------------------------------------------------------------------
        //Properties
        member this.dataDict with get()=DataDict
        member this.dataList with get()=DataList
        member this.groupDict with get()=GroupDict
        //--------------------------------------------------------------------
        //Local Methods
        //--------------------------------------------------------------------
        member this.SetDataValueByNameLocal (tag:string) (value:float) =
                lock(lockObj) (fun()->
                    try 
                        let i=DataDict.Item(tag)
                        DataList.Item(i).Value<-value
                        0
                    with
                    |_ -> 1    
                )
        //--------------------------------------------------------------------
        member this.SetDataStatusByNameLocal (tag:string) (status:int) =
                lock(lockObj) (fun()->
                    try 
                        let i=DataDict.Item(tag)
                        DataList.Item(i).Status<-status
                        0
                    with
                    |_ -> 1    
                )
        //--------------------------------------------------------------------
        member this.GetDataValueByNameLocal tag  = 
             lock(lockObj) (fun()->
                    try 
                        let i=DataDict.Item(tag)
                        DataList.Item(i).Value,0
                    with
                    |_ -> 0.0,1    
                )
        //--------------------------------------------------------------------
        //Group Methods
        member this.AddGroupLocal (group:string) (tags:string[])=
            if GroupDict.ContainsKey(group) then
                (1,[||])
            else
                let l= Array.zeroCreate(tags.Length)
                let e= Array.zeroCreate(tags.Length) //mark tags that doesn't exit
                //Check if all tags exit
                let mutable k=0
                let mutable j=0
                for i in tags do
                    if DataDict.ContainsKey(i) then
                        l.[k]<-DataDict.Item(i)
                    else
                        e.[k]<-1
                        j<-1
                    k<-k+1
                if j=0 then
                    //
                    GroupDict.TryAdd(group,l)|>ignore
                    (j,e)
                else
                    (j,e)
                
        //--------------------------------------------------------------------
        member this.DeleteGroupLocal (group:string)=
            if GroupDict.ContainsKey(group) then
                GroupDict.TryRemove(group)|>ignore
                0
            else
                1
        //--------------------------------------------------------------------
        member this.GetGroupLocal (group:string)=
            if GroupDict.ContainsKey(group) then
                let l=GroupDict.Item(group)
                let a= Array.zeroCreate(l.Length)
                let mutable k=0
                for i in l do
                    a.[k]<-DataList.[i]
                    k<-k+1
                (0,a)
            else
                (1,[||])
        //--------------------------------------------------------------------
        member this.GetGroupArrayLocal (group:string)=
            if GroupDict.ContainsKey(group) then
                let l=GroupDict.Item(group)
                let s= Array.zeroCreate(l.Length)
                let v : float array= Array.zeroCreate(l.Length)
                let mutable k=0
                for i in l do
                    v.[k]<-DataList.[i].Value                    
                    s.[k]<-DataList.[i].Status
                    k<-k+1
                (0,v,s)
            else
                (1,[||],[||])
        //--------------------------------------------------------------------    
        member this.RegisterDataLocal (tag:DataItem) = 
            lock(lockObj) (fun()->
                if DataDict.TryAdd(tag.Tag,DataList.Count) then
                    DataList.Add(tag)
                    0
                else
                    1)
        //--------------------------------------------------------------------
        member this.RegisterDataArrayLocal (tags:DataItem[])  = 
                //this.RegisterTag(i) 
            lock(lockObj) (fun()->
                if List.fold (&&) false [for t in tags ->DataDict.ContainsKey(t.Tag)] then
                    1
                else
                    for i in tags do
                        DataDict.TryAdd(i.Tag,DataList.Count)|>ignore
                        DataList.Add(i)
                    0)                
        //--------------------------------------------------------------------

        interface IDataServer with

            member this.AddGroup (group:string) (tags:string[])=
                   this.AddGroupLocal group tags
            //----------------------------------------------------------------
            member this.DeleteGroup group=
                   this.DeleteGroupLocal group
           //-----------------------------------------------------------------
            member this.GetGroup group=
                   this.GetGroupLocal group
            //----------------------------------------------------------------
            member this.SetDataValueByName tag value =
                this.SetDataValueByNameLocal tag value
            //----------------------------------------------------------------        
            member this.GetDataValueByName tag  =
                this.GetDataValueByNameLocal tag 
            //----------------------------------------------------------------
            member this.SetDataStatusByName tag status =
                this.SetDataStatusByNameLocal tag status
            //----------------------------------------------------------------
            member this.RegisterData (tag:DataItem) = 
                this.RegisterDataLocal(tag)
            //----------------------------------------------------------------
            member this.RegisterDataArray (tags:DataItem[])=
                this.RegisterDataArrayLocal(tags)
            //----------------------------------------------------------------
            member this.GetGroupArray (group:string)=
                this.GetGroupArrayLocal group
            //----------------------------------------------------------------

The file DS_Server.fsx implements the Data Server host.
DS_Server.fsx code:
First, the service is started at the server F# interactive window.

//Server Host
#r "System.ServiceModel"
#r "System.Runtime.Serialization"

//namespace RTE.DataServer

open System
open System.ServiceModel
open System.Collections.Concurrent
open System.Collections.Generic
open System.Runtime.Serialization

#I @"C:\Users\caxelrud\Documents\Visual Studio 2010\Projects\MailboxWCF"
#load "DS_Contract.fsx"

let DS= new DS_Contract.DataServer()


let h=
    //System.ServiceModel
    let baseAddress= new Uri("net.pipe://RTE/dataserver/ds1")
    let s = new ServiceHost(DS,[|baseAddress|])
    //Programmatic Binding
    let binding=
        let b=
            new NetNamedPipeBinding(securityMode=NetNamedPipeSecurityMode.None)
        b 
    s.AddServiceEndpoint(typeof,binding,baseAddress) |>ignore
    s

h.Open()

Now, we can use local methods to register and manipulate data items in the server F# interactive window.


// Test Local 
//Item
let tag1=new DS_Contract.DataItem(Tag="Tag_1",Description="my tag 1",Value=5.0);;
> val tag1 : DS_Contract.DataItem =
Tag=Tag_1, Description=my tag 1 Value=5.000000, Status=0

let r1=DS.RegisterDataLocal(tag1)
printfn "Tag_1 value=%A" (DS.GetDataValueByNameLocal("Tag_1"));;
> Tag_1 value=(5.0, 0)

let r2=DS.SetDataValueByNameLocal "Tag_1" 10.0
printfn "Tag_1 value=%A" (DS.GetDataValueByNameLocal("Tag_1"));;
> Tag_1 value=(10.0, 0)

printfn "tag 1=%f" (DS.dataList.Item(DS.dataDict.Item("Tag_1")).Value) //only local operation
> tag 1=10.000000

//Array
let tag2=new DS_Contract.DataItem(Tag="Tag_2",Description="my tag 2",Value=2.0)
let tag3=new DS_Contract.DataItem(Tag="Tag_3",Description="my tag 3",Value=3.0)
let tag4=new DS_Contract.DataItem(Tag="Tag_4",Description="my tag 4",Value=4.0)
let r3=DS.RegisterDataArrayLocal([|tag2;tag3;tag4|])
printfn "value=%A" (DS.GetDataValueByNameLocal("Tag_2")) 
> value=(2.0, 0)

//only local FSI
DS.dataDict;;
> val it : ConcurrentDictionary =
dict [("Tag_2", 1); ("Tag_4", 3); ("Tag_1", 0); ("Tag_3", 2)]

DS.dataList;;
> val it : List =
seq
[Tag=Tag_1, Description=my tag 1 Value=10.000000, Status=0
{Description = "my tag 1";
Status = 0;
Tag = "Tag_1";
Value = 10.0;};
Tag=Tag_2, Description=my tag 2 Value=2.000000, Status=0
{Description = "my tag 2";
Status = 0;
Tag = "Tag_2";
Value = 2.0;};
Tag=Tag_3, Description=my tag 3 Value=3.000000, Status=0
{Description = "my tag 3";
Status = 0;
Tag = "Tag_3";
Value = 3.0;};
Tag=Tag_4, Description=my tag 4 Value=4.000000, Status=0
{Description = "my tag 4";
Status = 0;
Tag = "Tag_4";
Value = 4.0;}]

//Group
let r4=DS.AddGroupLocal "G1" [|"Tag_1";"Tag_2"|];;
> val r4 : int * int [] = (0, [|0; 0|])

let r5=DS.AddGroupLocal "G2" [|"Tag_1";"Tag_2";"Tag_4"|]
> val r5 : int * int [] = (0, [|0; 0; 0|])

let r6=DS.GetGroupLocal "G1";;
> val r6 : int * DS_Contract.DataItem [] =
(0,
[|Tag=Tag_1, Description=my tag 1 Value=10.000000, Status=0;
Tag=Tag_2, Description=my tag 2 Value=2.000000, Status=0|])

printfn "transaction status=%d, value=%A" (fst r6) (snd r6);;
> transaction status=0, value=[|Tag=Tag_1, Description=my tag 1 Value=10.000000, Status=0;
Tag=Tag_2, Description=my tag 2 Value=2.000000, Status=0|]

let r7,r8,r9=DS.GetGroupArrayLocal "G1"
printfn "transaction status=%d, value=%A, status=%A" r7 r8 r9;;
> transaction status=0, value=[|10.0; 2.0|], status=[|0; 0|]

//only local FSI
DS.groupDict;;
> val it : ConcurrentDictionary =
dict [("G2", [|0; 1; 3|]); ("G1", [|0; 1|])]


//-------------------------------------------------------------------------
Console.WriteLine("Server is running. Press return to exit");;
Console.ReadLine()|>ignore
h.Close() //close it before updates !

DS_Client.fsx code:
The client is started at the client F# interactive window.

//Client
#r "System.ServiceModel"
#r "System.Runtime.Serialization"

open System
open System.ServiceModel
open System.Collections.Concurrent
open System.Collections.Generic
open System.Runtime.Serialization

#I @"C:\Users\caxelrud\Documents\Visual Studio 2010\Projects\MailboxWCF"
#load "DS_Contract.fsx"

let binding=new NetNamedPipeBinding(securityMode=NetNamedPipeSecurityMode.None)
let address=new EndpointAddress("net.pipe://RTE/dataserver/ds1")
let factory= new ChannelFactory(binding,address)
let channel=factory.CreateChannel()

// Remote Test
let tag10=new DS_Contract.DataItem(Tag="Tag_10",Description="my tag 10",Value=5.0);;
> val tag10 : DS_Contract.DataItem =
  Tag=Tag_10, Description=my tag 10 Value=5.000000, Status=0

let r1=channel.RegisterData(tag10)
printfn "Tag_10 value=%A" (channel.GetDataValueByName("Tag_10"));;
> Tag_10 value=(5.0, 0)

let r2=channel.SetDataValueByName "Tag_10" 10.0
printfn "Tag_10 value=%A" (channel.GetDataValueByName("Tag_10"));;
> Tag_10 value=(10.0, 0)

//Array
let tag20=new DS_Contract.DataItem(Tag="Tag_20",Description="my tag 20",Value=2.0)
let tag30=new DS_Contract.DataItem(Tag="Tag_30",Description="my tag 30",Value=3.0)
let tag40=new DS_Contract.DataItem(Tag="Tag_40",Description="my tag 40",Value=4.0)
let r3=channel.RegisterDataArray([|tag20;tag30;tag40|])
printfn "value=%A" (channel.GetDataValueByName("Tag_20"));; 
> value=(2.0, 0)

//Group
let r4=channel.AddGroup "G10" [|"Tag_10";"Tag_20"|];;
let r5=channel.AddGroup "G20" [|"Tag_10";"Tag_20";"Tag_40"|];;
let r6=channel.GetGroup "G10";;
printfn "transaction status=%d, value=%A" (fst r6) (snd r6);;
> transaction status=0, value=[|Tag=Tag_10, Description=my tag 10 Value=10.000000, Status=0;
  Tag=Tag_20, Description=my tag 20 Value=2.000000, Status=0|]

let r7,r8,r9=channel.GetGroupArray "G10"
printfn "transaction status=%d, value=%A, status=%A" r7 r8 r9;;
> transaction status=0, value=[|10.0; 2.0|], status=[|0; 0|]

Wednesday, December 15, 2010

Agent based example using F# MailboxProcessor

Introduction

The following simple example shows how to use F# Mailbox to develop an Agent (or Actor) based application.
Each Agent is a state-machine.
Each Agent do actions based on receiving messages as changing state (state-machine) or replying to the caller.
The action can be done asynchronous, allowing the application to easy scale keeping performance.

You can test and debug the application using F# interactive by send/receive messages to the Agents interactive.

Example

The example uses 3 Agents.

The agent2 and agent3 are of type agentCalc. It receives messages as 'Two' and 'Three' together with an integer. These Agents multiply the receiving number by 2 or by 3 depending on the message.

The agent1 has 2 states: idle and running. It also keeps a totalizer.
It receives messages as: 'Start', 'Exit', 'Reset', 'GoIdle', 'Next', 'Show' and 'State'.
These messages may have different meaning depending of the state of agent1.
For example:
In 'idle':
   - 'Start': change the state to 'running'
   - 'Exit': agent stop execution
In 'running':
    - 'GoIdle': change the state to 'idle'
    - 'Reset': totalizer goes to 0
    - 'Next': generates 2 random integers (0..1) and (0..9).
                The first number is used to decide if the message to be sent will be 'Two' or 'Three'.
                Then it send a message to agent2 and agent3 passing the second number.
                The returning number from agent2 and agent3 is added to the totalizer.

The following script implements the agents:
open System

type msg1 =
    | Start
    | Exit
    | Reset
    | GoIdle
    | Next
    | Show of AsyncReplyChannel
    | State of AsyncReplyChannel


type msg2 =
    Two of int*AsyncReplyChannel | Three of int*AsyncReplyChannel
    


let agentCalc ()=
    MailboxProcessor.Start(fun i ->
        let rec running =
            async {
                let! msg= i.Receive()
                match msg with
                | Two (n,replyChannel) ->
                    do replyChannel.Reply(n*2)
                    return! running
                | Three (n, replyChannel) ->
                    do replyChannel.Reply(n*3)
                    return! running
                }
        running )

let agent2= agentCalc ()
let agent3= agentCalc ()

let mutable total=0
let agent1=
        MailboxProcessor.Start(fun i ->
            let rec idle (total:int)=
                async {
                    let! msg= i.Receive()
                    match msg with
                    | Reset ->
                        let total=0
                        return! idle total
                    | Exit ->
                        return ()
                    | Start ->
                        return! running total
                    | Show (replyChannel)->
                        do replyChannel.Reply(total)
                        return! idle total
                    | State (replyChannel)->
                        do replyChannel.Reply("idle")
                        return! idle total
                    | _->
                        return! idle total                        
                    }
            and running total= 
                async{
                    let! msg= i.Receive()
                    match msg with
                    | Next ->
                        let z=new Random()
                        let n1=z.Next(1)
                        let n2=z.Next(9)
                        let total=
                                if n1=0 then
                                    let r1=agent2.PostAndReply(fun replyChannel->Two(n2,replyChannel))
                                    let r2=agent3.PostAndReply(fun replyChannel->Two(n2,replyChannel))
                                    total+r1+r2
                                else
                                    let r1=agent2.PostAndReply(fun replyChannel->Three(n2,replyChannel))
                                    let r2=agent3.PostAndReply(fun replyChannel->Three(n2,replyChannel))
                                    total+r1+r2
                        return! running total
                    | GoIdle ->
                        return! idle total
                    | Show (replyChannel) ->
                        do replyChannel.Reply(total)
                        return! running total
                    | State (replyChannel)->
                        do replyChannel.Reply("running")
                        return! running total
                    | Reset ->
                        let total=0
                        return! running total
                    | _->
                        return! idle total                        
                }
            idle total)


The following script tests interactivelly the agents:

let r1 = agent1.PostAndReply(fun reply -> State(reply))
printfn "State %s" r1 |>ignore;;
State idle 

let r2 = agent1.PostAndReply(fun reply -> Show(reply))
printfn "Total %d" r2 |>ignore;;
Total 0 

agent1.Post(Start)
let r3 = agent1.PostAndReply(fun reply -> State(reply))
printfn "State %s" r3 |>ignore;;
State running 

for i in 1..5 do
    agent1.Post(Next)
    let r4 = agent1.PostAndReply(fun reply -> Show(reply))
    printfn "Total %d" r4 |>ignore;;

Total 24
Total 48
Total 72
Total 96
Total 120

agent1.Post(Reset)
let r5 = agent1.PostAndReply(fun reply -> Show(reply))
printfn "Total %d" r5 |>ignore;;
Total 0 

agent1.Post(GoIdle)
let r6 = agent1.PostAndReply(fun reply -> State(reply))
printfn "State %s" r6 |>ignore;;
State idle 

agent1.Post(Exit);;

Monday, September 27, 2010

F# parallel programming using Agents

Consider a parallel program as a series of independent components that send and receive messages.This is referred as the agent model.

When an actor receives a message, the message is placed in a queue until the agent is ready to process it. When an agent process a message, it makes a decision about what to do with it based on its internal state and contents of the message. It might do processing, send a reply, send a new message for a different agent, etc.

F# provides the MailboxProcessor class for agent processing using Asynchronous Programming or more specific, Asynchronous Workflows.

Asynchronous Programming describes programs and operations that once started are executed in background and terminate at some later time.

Asynchronous Workflow allow you to perform asynchronous operations without the need for explicit callbacks. So you can write the code as if it were using synchronous execution, but actually, the code will execute asynchronously, supending and resuming as asynchronous operation complete.

In Asynchronous Programming, you want to avoid blocking of threads. Threads are expensive because allocate a 1 MB stack and other operational system implications.
In performance critical code, it is important to keep the number of blocked threads low. Ideally, you will only have as many as you have processors and you will have no blocked threads.

There are several topics related to the usage of MailboxProcessor. In this document I will analyze the handling of agent's persisted variables.

The code for the first agent example is:

let Act_RandomAccum= 
    MailboxProcessor.Start(fun (i:MailboxProcessor) ->
        let rec loop state =
            async {
                    let! r= i.Receive()
                    let z=new Random()
                    let v=z.NextDouble()
                    printfn "state=%f" (state+v)
                    return! loop(state+v)}
        loop 0.0)
(*
state=0.761920
state=1.613361
state=1.717945
state=1.912050
...
*)

Let's create an function that will run continuously and it will post messages to agents every 3 seconds (asynchronously). You can cancel it by executing Async.CancelDefaultToken().

let Act_clk=
    async {
            while true do
                do Act_RandomAccum.Post(1)
                //do Act_3States.Post(1)
                //do Act_3States_2.Post(1)
                do! Async.Sleep(3*1000)
           }

let cancelHandler(ex : OperationCanceledException)=
    printfn "Task canceled"

Async.TryCancelled(Act_clk, cancelHandler) 
|> Async.Start

Notice how variables are keep by being a input to the loop.
This is not covinient when you want to preserve several states.

The next snippet uses a record to preserve several variables. Change the comment lines in Act_clk to send messages to the following agent.

type States= {mutable s1:float; mutable s2:float; mutable s3:float; mutable sum:float; mutable Status: string}

let States1={s1=0.0;s2=0.0;s3=0.0;sum=0.0;Status="good"}

let Act_Random3States= 
    MailboxProcessor.Start(fun (i:MailboxProcessor) ->
        let rec loop (s:States) =
            async {
                    let! r= i.Receive()
                    let z=new Random()
                    let v=z.NextDouble()
                    s.s3<-s.s2; s.s2<-s.s1
                    s.s1<-v
                    s.sum<-s.s1+s.s2+s.s3
                    printfn "states,sum=%f %f %f %f" s.s1 s.s2 s.s3 s.sum
                    return! loop(s)}
        loop States1)
(* 
states,sum=0.023424 0.000000 0.000000 0.023424
states,sum=0.046045 0.023424 0.000000 0.069468
states,sum=0.068666 0.046045 0.023424 0.138134
states,sum=0.568861 0.068666 0.046045 0.683572
states,sum=0.591482 0.568861 0.068666 1.229009
states,sum=0.614103 0.591482 0.568861 1.774447
states,sum=0.636724 0.614103 0.591482 1.842310
...
*)

The last example is more powerfull. The states are preserved in a class. Objects (as the random one) can be stored and don't need to be recreated in every execution.
Instance and static function can be defined (as fShift and fSum) to be used by the agent.
In other words, the agent defition is encapsulated in a class and its execution is done in the MailboxProcessor.

type StatesClass(s1:float,s2:float,s3:float,status:string)=
    let mutable _s1=s1 
    let mutable _s2=s2 
    let mutable _s3=s3 
    let mutable _status=status 
    let r=new Random()
    let mutable sum= StatesClass.fSum _s1 _s2 _s3  

    member this.Sum with get()=sum and set x=sum<-x
    member this.S1 with get()= _s1 and set x=_s1<-x
    member this.S2 with get()= _s2 and set x=_s2<-x
    member this.S3 with get()= _s3 and set x=_s3<-x
    member this.Status with get()= _status and set x=_status<-x
    member this.rndObj=r    
    static member fSum x y z=x+y+z
    member this.fShift()=
        this.S3<-this.S2; this.S2<-this.S1

let sc1=StatesClass(s1=0.0,s2=0.0,s3=0.0,status="ok")

let Act_3States_2= 
    MailboxProcessor.Start(fun (i:MailboxProcessor) ->
        let rec loop (s:StatesClass) =
            async {
                    let! r= i.Receive()
                    s.fShift()
                    let v=s.rndObj.NextDouble()
                    s.S1<-v
                    s.Sum <- StatesClass.fSum s.S1 s.S2 s.S3
                    printfn "states,sum=%f %f %f %f" s.S1 s.S2 s.S3 s.Sum
                    return! loop(s)}
        loop sc1)
(*
states,sum=0.396605 0.000000 0.000000 0.396605
states,sum=0.029249 0.396605 0.000000 0.425854
states,sum=0.197748 0.029249 0.396605 0.623602
states,sum=0.863481 0.197748 0.029249 1.090479
states,sum=0.641731 0.863481 0.197748 1.702961
states,sum=0.867453 0.641731 0.863481 2.372665
states,sum=0.880071 0.867453 0.641731 2.389254
...
*)

Thursday, September 23, 2010

Large-scale nonlinear optimization for F#

This blog shows how to use Ipopt (https://projects.coin-or.org/Ipopt) in F#.
Ipopt (Interior Point OPTimizer) is a software package for large-scale nonlinear optimization.
It is designed to find (local) solutions of mathematical optimization problems of the from
   min f(x)
x in R^n
s.t.    g_L <= g(x) <= g_U
        x_L <= x <= x_U

where f(x): R^n --> R is the objective function, and g(x): R^n --> R^m are the constraint functions. The vectors g_L and g_U denote the lower and upper bounds on the constraints, and the vectors x_L and x_U are the bounds on the variables x. The functions f(x) and g(x) can be nonlinear and nonconvex, but should be twice continuously differentiable. Note that equality constraints can be formulated in the above formulation by setting the corresponding components of g_L and g_U to the same value.

The binaries can be downloaded from http://www.coin-or.org/download/binary/Ipopt/. I am using the library Ipopt38.dll.
The interface for .NET can be downloaded from http://code.google.com/p/csipopt/.
Following the instruction on README.txt, you can create a class library Cureos.Numeric.dll by running the batch script build_library.bat.
Based on the C# example, provided by csipopt, I created the following F# code. I created a type, following the example, but it is not really necessary.

#I @"C:\Users\caxelrud\Documents\Visual Studio 2010\Projects\hs071_fsx"
#r "Cureos.Numerics.dll"
//Need Ipopt38.dll in path

open System
open System.IO
open Cureos.Numerics
open System.Runtime.InteropServices

type HS071=
    val public n : int
    val public m : int
    val public nele_jac : int
    val public nele_hess : int
    val public x_L : double[]
    val public x_U : double[]
    val public g_L : double[]
    val public g_U : double[]
    //set the number of variables and allocate space for the bounds
    //set the values for the variable bounds using n, x_L, x_U
    //Set the number of constraints and allocate space for the bounds using m, g_L, g_U
    //Set the number of nonzeros in the Jacobian of the constraints using nele_jac
    //Set the number of nonzeros in the Hessian of the Lagrangian (lower or upper triangual part only) using nele_hess

    new ()= {n=4;m=2;nele_jac = 8;nele_hess = 10;x_L =[1.0; 1.0; 1.0; 1.0 ];x_U = [5.0; 5.0; 5.0; 5.0];g_L=[25.0; 40.0];g_U = [Ipopt.PositiveInfinity; 40.0]}
 
    member public this.eval_f (n:int, x:double[], new_x:bool, [] obj_value:byref)=
        obj_value <- x.[0]*x.[3]*(x.[0]+x.[1]+x.[2])+x.[2]
        true

    member public this.eval_grad_f (n:int, x:double[], new_x:bool, grad_f:double[])=
        grad_f.[0] <- x.[0] * x.[3] + x.[3] * (x.[0] + x.[1] + x.[2])
        grad_f.[1] <- x.[0] * x.[3]
        grad_f.[2] <- x.[0] * x.[3] + 1.0
        grad_f.[3] <- x.[0] * (x.[0] + x.[1] + x.[2])
        true

    member public this.eval_g (n:int, x:double[], new_x:bool, m:int, g:double[])=

        g.[0] <- x.[0] * x.[1] * x.[2] * x.[3];
        g.[1] <- x.[0] * x.[0] + x.[1] * x.[1] + x.[2] * x.[2] + x.[3] * x.[3];
        true;

    member public this.eval_jac_g (n:int, x:double[], new_x:bool, m:int, nele_jac:int, iRow:int[], jCol:int[], values:double[])=
        if values = null then
            // set the structure of the jacobian
           // this particular jacobian is dense
            iRow.[0] <- 0
            jCol.[0] <- 0
            iRow.[1] <- 0
            jCol.[1] <- 1
            iRow.[2] <- 0
            jCol.[2] <- 2
            iRow.[3] <- 0
            jCol.[3] <- 3
            iRow.[4] <- 1
            jCol.[4] <- 0
            iRow.[5] <- 1
           jCol.[5] <- 1
           iRow.[6] <- 1
           jCol.[6] <- 2
           iRow.[7] <- 1
           jCol.[7] <- 3
       else
          // return the values of the jacobian of the constraints
          values.[0] <- x.[1] * x.[2] * x.[3] // 0,0
          values.[1] <- x.[0] * x.[2] * x.[3] // 0,1
          values.[2] <- x.[0] * x.[1] * x.[3] // 0,2
          values.[3] <- x.[0] * x.[1] * x.[2] // 0,3
          values.[4] <- 2.0 * x.[0] // 1,0
          values.[5] <- 2.0 * x.[1] // 1,1
          values.[6] <- 2.0 * x.[2] // 1,2
          values.[7] <- 2.0 * x.[3] // 1,3
       true

    member public this.eval_h (n:int, x:double[], new_x:bool, obj_factor:double,m:int, lambda:double[],  new_lambda:bool,nele_hess:int, iRow:int[], jCol:int[],values:double[]) =

        if values = null then
            // set the Hessian structure. This is a symmetric matrix, fill the lower left
            // triangle only.
            // the hessian for this problem is actually dense
            let mutable idx= 0; // nonzero element counter
            for row in 0..3 do
                 for col in 0..row do
                      iRow.[idx] <- row
                      jCol.[idx] <- col
                      idx<-idx+1
        else
            // return the values. This is a symmetric matrix, fill the lower left
            // triangle only
            // fill the objective portion
            values.[0] <- obj_factor * (2.0 * x.[3]); // 0,0
            values.[1] <- obj_factor * (x.[3]); // 1,0
            values.[2] <- 0.0; // 1,1
            values.[3] <- obj_factor * (x.[3]); // 2,0
            values.[4] <- 0.0; // 2,1
            values.[5] <- 0.0; // 2,2
            values.[6] <- obj_factor * (2.0 * x.[0] + x.[1] + x.[2]); // 3,0
            values.[7] <- obj_factor * (x.[0]); // 3,1
            values.[8] <- obj_factor * (x.[0]); // 3,2
            values.[9] <- 0.0; // 3,3
            // add the portion for the first constraint
            values.[1] <- values.[1]+ lambda.[0] * (x.[2] * x.[3]); // 1,0 */
            values.[3] <- values.[3]+ lambda.[0] * (x.[1] * x.[3]); // 2,0 */
            values.[4] <- values.[4]+ lambda.[0] * (x.[0] * x.[3]); // 2,1 */
            values.[6] <- values.[6]+ lambda.[0] * (x.[1] * x.[2]); // 3,0 */
            values.[7] <- values.[7]+ lambda.[0] * (x.[0] * x.[2]); // 3,1 */
            values.[8] <- values.[8]+ lambda.[0] * (x.[0] * x.[1]); // 3,2 */
            // add the portion for the second constraint
            values.[0] <- values.[0]+ lambda.[1] * 2.0; // 0,0 */
            values.[2] <- values.[2]+ lambda.[1] * 2.0; // 1,1 */
            values.[5] <- values.[5]+ lambda.[1] * 2.0; // 2,2 */
            values.[9] <- values.[9]+ lambda.[1] * 2.0; // 3,3 */
        true

// create the IpoptProblem

let p = new HS071()
// allocate space for the initial point and set the values
let x = [1.0; 5.0; 5.0; 1.0 ]
let mutable obj=0.0
let d_eval_g= new EvaluateConstraintsDelegate(fun (n:int) (x:double[]) (new_x:bool) (m:int) (g:double[])-> p.eval_g(n,x,new_x,m,g))
let d_eval_f= new EvaluateObjectiveDelegate(fun n x new_x obj_value-> p.eval_f(n,x,new_x,&obj_value))
let d_eval_grad_f= new EvaluateObjectiveGradientDelegate(fun n x new_x grad_f-> p.eval_grad_f(n,x,new_x,grad_f))
let d_eval_jac_g= new EvaluateJacobianDelegate(fun n x new_x m nele_jac iRow jCol values-> p.eval_jac_g(n,x,new_x,m,nele_jac,iRow,jCol,values))
let d_eval_h= new EvaluateHessianDelegate(fun n x new_x obj_factor m lambda new_Lambda nele_hess iRow jCol values-> p.eval_h(n,x,new_x,obj_factor,m,lambda,new_Lambda,nele_hess,iRow,jCol,values))

let problem = new Ipopt(p.n, p.x_L, p.x_U, p.m, p.g_L, p.g_U, p.nele_jac, p.nele_hess, d_eval_f, d_eval_g, d_eval_grad_f, d_eval_jac_g, d_eval_h)

// Set some options. Note the following ones are only examples, they might not be suitable for your problem.
problem.AddOption("tol", 1e-7) >ignore
problem.AddOption("mu_strategy", "adaptive") >ignore
problem.AddOption("output_file", "hs071_fs.txt") >ignore

// solve the problem
let status=problem.SolveProblem(x, &obj, null, null, null, null);
printfn "Optimization return status:%A" status
    for i in 0..3 do
        printfn "x[%d]=%f" i x.[i]
(*
Optimization return status:Solve_Succeeded
x[0]=1.000000
x[1]=4.743000
x[2]=3.821150
x[3]=1.379408
*)

Tuesday, September 14, 2010

Graphics in F# interactive

The following script snippets present graphics functionalities available in FlyingFrog.Graphics namespace from FSharpForVisualization assembly.
This functionalities allow the usage of F# interactive with plots much like other tools as Matlab or Python with MatPlotLib.
#I @"C:\cs\Docs"
#r "FSharp.PowerPack.dll"

// Reference WPF
#r "PresentationCore.dll"
#r "PresentationFramework.dll"
#r "WindowsBase.dll"

// Reference to FlyingFrog
#r @"FSharpForVisualization.dll"


open System.Windows
open System.Windows.Media
open FlyingFrog.Graphics;;

//Plot
Plot([Function sin], (-6., 6.));;

// Plot coordinates
let xys =
   [ for x in -5 .. 5 ->
       let x = float x
       x, sin x ]

Plot([Function sin; Data xys], (-6., 6.), Labels=(Text "x", Text "sin(x)"));;


// Progressive creation of a plot. Begin with an empty plot:
let g = Plot([], (-6., 6.));;

// And then specify a function to be plotted:
g.Data <- [Function sinc; Data[System.Math.PI, 0.]];;

// Specify as many functions:
g.Data <- [Function sin; Function cos];;

// Alter the styles:
g.Styles <- [Brushes.Black; Brushes.Gray];;



// 3D surface plot
let f x z =
    let r = 5. * sqrt(x*x + z*z)
    sin(r + 3. * x) / r

Plot3D(f, (-5., 5.), (-5., 5.), (-0.1, 1.));;
 
 
// Vector graphic shapes are styled by a sequence of strokes or fills:

let styles =

    [ Stroke(Brushes.Red, {width = 0.001}) ]
let geometry =

    let n = 50
    [ for i in 0 .. n do
         for j in 0 .. 1 do
            let x, y = float i / float n, float j
            yield Contour.Open(Point(x, y), [ Segment.Line(Point(y, 1. - x)) ]) ]
let scene = Shape(styles, geometry)

View scene;;



Plot([], (0., 1.), (0., 1.), Epilog=scene);;
 


// 3D object - icosahedron
let scene1 = Shape3D(Polyhedra.Icosahedron.mesh, Brushes.White)
View3D scene1;;

Sunday, September 12, 2010

Vector, Matrix and Linear Algebra in F#

The following script snippets present the vector and matrix functionalities in Microsoft.Fsharp.Numerics namespace from Fsharp.PowerPack assemble.
Also, they present the Linear Algebra functionalities from FlyingFrog namespace from FSharpForNumeric assemble.

//F# documentation at
//http://msdn.microsoft.com/library/dd233154(VS.100).aspx
//Fsharp.PowerPack documentation at
//http://research.microsoft.com/en-us/um/cambridge/projects/fsharp/manual/namespaces.html
//Flying frogs numerics documentation at
//http://www.ffconsultancy.com/products/fsharp_for_numerics/docs/namespaces.html


#I @"C:\Program Files (x86)\FSharpPowerPack-2.0.0.0\bin"
#r "Fsharp.PowerPack.dll"
#r "FSharpForNumerics.dll"
#r "FSharpForVisualization.dll"

open Microsoft.FSharp.Math
open FlyingFrog
open FlyingFrog.Graphics
open FlyingFrog.Graphics.Typeset;;

//vector (float)
let v1 = vector [1.0;1.0;1.0]
let v2 = vector [2.0;2.0;2.0]
let v3 = rowvec [3.0;3.0;3.0];;
val v1 : vector = vector [|1.0; 1.0; 1.0|]
val v2 : vector = vector [|2.0; 2.0; 2.0|]
val v3 : rowvec = rowvec [|3.0; 3.0; 3.0|]

//Simple operations
2.0*v1;;
v1+v2;;
v1-v2;;
v1.*v2;; //per element multiplication
let m1 = v2*v3;;

val it : Vector = vector [|2.0; 2.0; 2.0|]
val it : Vector = vector [|3.0; 3.0; 3.0|]
val it : Vector = vector [|-1.0; -1.0; -1.0|]
val it : Vector = vector [|2.0; 2.0; 2.0|]
val m1 : Matrix = matrix [[6.0; 6.0; 6.0]
[6.0; 6.0; 6.0]
[6.0; 6.0; 6.0]]
//Transpose
v1.Transpose;;
val it : RowVector = rowvec [|1.0; 1.0; 1.0|]

//Vector Generic
type intVec = Vector
let intVec ll = Vector.Generic.ofSeq ll
let A = intVec [1;2;3;4]
let B = intVec [1;1;1;1]
A+B;;
A*10;;
Vector.Generic.sum A;;

type intVec = Vector
val intVec : seq<'a> -> Vector<'a>
val A : Vector = vector [|1; 2; 3; 4|]
val B : Vector = vector [|1; 1; 1; 1|]

>val it : Vector = vector [|2; 3; 4; 5|]
> val it : Vector = vector [|10; 20; 30; 40|]
> val it : int = 10

//matrix (float)
let A = matrix [ [ 1.0; 7.0; 2.0 ];
[ 1.0; 3.0; 1.0 ];
[ 2.0; 9.0; 1.0 ]; ]
let B = matrix [ [ 10.0; 70.0; 20.0 ];
[ 10.0; 30.0; 10.0 ];
[ 20.0; 90.0; 10.0 ]; ];;

val A : matrix = matrix [[1.0; 7.0; 2.0]
[1.0; 3.0; 1.0]
[2.0; 9.0; 1.0]]
val B : matrix = matrix [[10.0; 70.0; 20.0]
[10.0; 30.0; 10.0]
[20.0; 90.0; 10.0]]
//print
let f (i:Matrix<_>)= printfn "%A" (i.ToArray2D())
let g (i:Vector<_>)= printfn "%A" i;;

val f : Matrix<'a> -> unit
val g : Vector<'a> -> unit

// sum, subtract
f (A+B)
> [[11.0; 77.0; 22.0]
[11.0; 33.0; 11.0]
[22.0; 99.0; 11.0]]
> [[-9.0; -63.0; -18.0]
[-9.0; -27.0; -9.0]
[-18.0; -81.0; -9.0]]

// matrix product
f (A*B);;

[[120.0; 460.0; 110.0]
[60.0; 250.0; 60.0]
[130.0; 500.0; 140.0]]

// element-wise product
f (A.*B);;

[[10.0; 490.0; 40.0]
[10.0; 90.0; 10.0]
[40.0; 810.0; 10.0]]
// scalar product
f (A * 2.0);;

[[2.0; 14.0; 4.0]
[2.0; 6.0; 2.0]
[4.0; 18.0; 2.0]]
f (2.0 * A);;
[[2.0; 14.0; 4.0]
[2.0; 6.0; 2.0]
[4.0; 18.0; 2.0]]
// negation of a matrix
f (-A);;

[[-1.0; -7.0; -2.0]
[-1.0; -3.0; -1.0]
[-2.0; -9.0; -1.0]]

// matrix-vector product
let b = vector [5.;8.;9.];;
g (A*b);;

val b : vector = vector [|5.0; 8.0; 9.0|]
> vector [|79.0; 38.0; 91.0|]

// Dimension
let dim = A.Dimensions;;

val dim : int * int = (3, 3)

//Transpose
Let A' = A.Transpose;;

val A' : Matrix = matrix [[1.0; 1.0; 2.0]
[7.0; 3.0; 9.0]
[2.0; 1.0; 1.0]]
//Number of rows,columns
let nrow = A.NumRows
let ncol = A.NumColslet;;

val nrow : int = 3
val ncol : int = 3

//Copy
let Anew = A.Copy();;

val Anew : Matrix = matrix [[1.0; 7.0; 2.0]
[1.0; 3.0; 1.0]
[2.0; 9.0; 1.0]]
// convert to a Array2D type
let Aarr = A.ToArray2D();;

val Aarr : float [,] = [[1.0; 7.0; 2.0]
[1.0; 3.0; 1.0]
[2.0; 9.0; 1.0]]
//Convert to vector (when the matrix has one row or column)
let Avec = A.ToVector()
let Arvec = A.ToRowVector();;

val Avec : Vector = vector [|1.0; 1.0; 2.0|]
val Arvec : RowVector = rowvec [|1.0; 7.0; 2.0|]

//operator [,]
A.[2,2];;
A.[2,2] <- 100.0;; A.[2,2];;

> val it : float = 100.0
> val it : unit = ()
> val it : float = 100.0

//sub-matrix
A.[1..2,1..2];;
A.[1..2,1..2] <- matrix [[2.;3.]; [8.;9.;]];; A.[1..2,1..2];;

> val it : Matrix = matrix [[3.0; 1.0]
[9.0; 1.0]]
> val it : unit = ()
> val it : Matrix = matrix [[2.0; 3.0]
[8.0; 9.0]]
//columns, rows
A.Column 2;;
A.Row 2;;
A.Columns (1,2);;
A.Rows (1,2);;

> val it : Vector = vector [|2.0; 3.0; 9.0|]
> val it : RowVector = rowvec [|2.0; 8.0; 9.0|]
> val it : Matrix = matrix [[7.0; 2.0]
[2.0; 3.0]
[8.0; 9.0]]
> val it : Matrix = matrix [[1.0; 2.0; 3.0]
[2.0; 8.0; 9.0]]
// sum and prod of all elements
let Asum = Matrix.sum A;;
let Aprod = Matrix.prod A;;

> val Asum : float = 35.0
> val Aprod : float = 12096.0

// create a matrix with 1
let C = Matrix.create 5 5 1.0;;

> val C : matrix = matrix [[1.0; 1.0; 1.0; 1.0; 1.0]
[1.0; 1.0; 1.0; 1.0; 1.0]
[1.0; 1.0; 1.0; 1.0; 1.0]
[1.0; 1.0; 1.0; 1.0; 1.0]
[1.0; 1.0; 1.0; 1.0; 1.0]]

// create a matrix with a function
let table = Matrix.init 3 3 (fun i j ->; (float i + 1.) * (float j + 1.));;

> val table : matrix = matrix [[1.0; 2.0; 3.0]
[2.0; 4.0; 6.0]
[3.0; 6.0; 9.0]]
// Identity
let I = Matrix.identity 2;;
val I : matrix = matrix [[1.0; 0.0]
[0.0; 1.0]]
// trace
let Atrace = Matrix.trace A;;

val Atrace : float = 12.0

//Map a function to the elements
let Asqr = Matrix.map (fun x -> x*x) A;;

> val Asqr : matrix = matrix [[1.0; 49.0; 4.0]
[1.0; 4.0; 9.0]
[4.0; 64.0; 81.0]]
//Sparse matrix
let D = Matrix.initSparse 3 3 [ (0,0,1.0); (1,1,2.0); (2,2,3.0); ];;

> val D : matrix = matrix [[1.0; 0.0; 0.0]
[0.0; 2.0; 0.0]
[0.0; 0.0; 3.0]]
D.IsSparse;;
> val it : bool = true

//Matrix Generic
type intmatrix = Matrix
let intmatrix ll = Matrix.Generic.ofSeq ll
let C = intmatrix [ [1;2]; [3;4] ]
let D = intmatrix [ [1;1]; [1;1] ]
Matrix.Generic.inplaceAdd C D;;
C+D;;
C * 10;;
type intmatrix = Matrix
val intmatrix : seq<#seq<'b>> -> Matrix<'b>
val C : Matrix = matrix [[2; 3]
[4; 5]]
val D : Matrix = matrix [[1; 1]
[1; 1]]
> val it : Matrix = matrix [[3; 4]
[5; 6]]
> val it : Matrix = matrix [[20; 30]
[40; 50]]
//Linear Algebra - float

//Solve linear equations
let A1 =
[ [1; 2; 3];
[1; 4; 9];
[1; 8; 27] ]
let A2 = Matrix.ofSeq(Seq.map (Seq.map float) A1)
let b = vector [2.; 3.; 4.]
let x = LinearAlgebra.Float.linear_solve A b;;

val A1 : int list list = [[1; 2; 3]; [1; 4; 9]; [1; 8; 27]]
val A2 : matrix = matrix [[1.0; 2.0; 3.0]
[1.0; 4.0; 9.0]
[1.0; 8.0; 27.0]]
val b : vector = vector [|2.0; 3.0; 4.0|]
val x : vector = vector [|0.5; 1.0; -0.1666666667|]

//Determinant
LinearAlgebra.Float.determinant A2;;

> val it : float = -12.0

//Inverse
LinearAlgebra.Float.inverse A2;;
> val it : matrix = matrix [[3.0; -2.5; 0.5]
[-1.5; 2.0; -0.5]
[0.3333333333; -0.5; 0.1666666667]]
//svd
LinearAlgebra.Float.svd A2;;

> val it : matrix * matrix * matrix =
(matrix [[0.1166955507; 0.6917676817; 0.7126286712]
[0.3269819117; 0.6507677484; -0.6852621157]
[0.9377979409; -0.3129837252; 0.1502538182]],
matrix [[30.0404326; 0.0; 0.0]
[0.0; 1.878075836; 0.0]
[0.0; 0.0; 0.2126972812]],
matrix [[0.0459872007; 0.5481949583; 0.835085304]
[0.301051; 0.7894977011; -0.5348473383]
[0.9524985421; -0.2759993977; 0.1287278512]])

//Linear Algebra – rational

//Solve linear equations
let A3 =
[ [1; 2; 3];
[1; 4; 9];
[1; 8; 27] ]
let A4 = Matrix.Generic.ofSeq(Seq.map (Seq.map BigNum.FromInt) A3)
let b1 = Vector.Generic.ofSeq [2N; 3N; 4N]
let x1 = LinearAlgebra.Rational.linear_solve A4 b1;;
A4 * x1;;

> val A3 : int list list = [[1; 2; 3]; [1; 4; 9]; [1; 8; 27]]
val A4 : Matrix = matrix [[1N; 2N; 3N]
[1N; 4N; 9N]
[1N; 8N; 27N]]
val b1 : Vector = vector [|2N; 3N; 4N|]
val x1 : Vector = vector [|1/2N; 1N; -1/6N|]
> val it : Vector = vector [|2N; 3N; 4N|]


//Inverse - rational
LinearAlgebra.Rational.inverse A4;;

> val it : Matrix = matrix [[3N; -5/2N; 1/2N]
[-3/2N; 2N; -1/2N]
[1/3N; -1/2N; 1/6N]]

//Linear Algebra - complex

type cmatrix = Matrix
let cmatrix ll = Matrix.Generic.ofSeq ll

let A5 = cmatrix [ [ complex 1.0 1.0; complex 1.0 1.0];
[ complex 1.0 1.0; complex 1.0 1.0]];;

> type cmatrix = Matrix
val cmatrix : seq<#seq<'b>> -> Matrix<'b>
val A5 : Matrix = matrix [[1r+1i; 1r+1i]
[1r+1i; 1r+1i]]
//Determinant - complex
LinearAlgebra.Complex.determinant A5;;

> val it : Complex = 0r+0i {Conjugate = 0r+0i;
ImaginaryPart = 0.0;
Magnitude = 0.0;
Phase = 0.0;
RealPart = 0.0;
i = 0.0;
r = 0.0;}
f (A-B);;

Saturday, September 11, 2010

F#- Introduction

F# is a first class language in .NET and Visual Studio together with C#, VB and C++. Solutions are created and debugged in the same way. IntelliSense is available. Generated assembles can be used in the .NET solutions. Execution speed is equivalent to C#.

But, F# also has unique features as F# Interactive and F# Scripting Files. These features are attractive to scientific and engineering developers that are used to tools like Matlab or Python.

Installation

A free installation is available as follow:

1) Install .NET 4.0 from
http://www.microsoft.com/downloads/en/details.aspx?FamilyID=9cfb2d51-5ff4-4491-b0e5-b386f32c0992&displaylang=en

2) Install Visual Studio Shell 10 from

http://www.microsoft.com/downloads/en/details.aspx?familyid=8E5AA7B6-8436-43F0-B778-00C3BCA733D3&displaylang=en

3) Install F# 2.0 from

http://www.microsoft.com/downloads/en/details.aspx?FamilyID=f8c623ae-aef6-4a06-a185-05f59be47d67&displaylang=en

4) Add extra content from:

http://fsharppowerpack.codeplex.com/releases/view/40168

Using F# Interactive

In order to open an interactive session, start the following executable:

C:\Program Files (x86)\Microsoft F#\v4.0\fsi.exe

Interactive session is also available inside Visual Studio. You can be shown it with ‘View/Other Windows/F# Interactive’.

Using F# Script Files

Script files have extension ‘.fsx’. You can create them with Notepad. By right-clicking on Windows Explorer, you can execute it by selecting “Run with F# Interactive”. This will start the script file with FSI.

Using Visual Studio to create script files has the advantage of IntelliSense. Create an script by using ‘File/New/File/F# Script File’.

Compiling F#

In order to generate an executable or a library, create a ‘.fs’ file in Notepad and use it with:

C:\Program Files\Microsoft F#\v4.0\fsc.exe

Using Visual Studio to create source files has the advantage of IntelliSense and easy reference to other assembles and .NET resources. Create a project by using ‘File/New/Project/F# Application or F# Library’.

Using Source File and Interactive

Code can be entered in Visual Studio source file. You can transfer and execute sections of the code by marking it and pressing ‘Alt+Enter’. This is a method to debug the code, traditionally used by systems with interactive capability.