Calling Matlab from a Scipt

From OpenSeesWiki
Revision as of 01:53, 12 March 2011 by Fmk (talk | contribs)
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to navigation Jump to search

It is possible to call matlab from a script to perform some calculation needed in the script. It is made possible by using the -r flag among others when invoking matlab from the command line (dos shell). The following is a little example demonstrating how this is done.


The example is a simple truss example where I am trying to determine the value of the parameter A needed to hit a target displacement.

  1. I start off with setting some variables, my target displacement, a tolerance, and two bounding values on my parameter of interest Amin and Amax.
  2. I then define a proc getA(). In this procedur I write the two bounding values to a file, I invoke a matlab script to determine a new A, and then I read the result from a file.
  3. Next I create my truss model and the analysis and perform a single analysis step
  4. Then I iterate until i have converged on a solution or hit a bounding value. In the iteration I call the matlab function to determine a new A, set the parameter in the elements, do an analysis step to obtain new displacement, and check versus target dispacement
  5. Finally I print a message indicating succcess or failure
#                                                                                                   
# set some variables                                                                                
#                                                                                                   

set targetDisp 1.2; # target displacement                                                           
set tol 0.01;       # tolerance                                                                     

set Amin 1.0
set Amax 10.0

#                                                                                                   
# define a procedure: getA() which calls matlab script getA.m in which our new A is computed           
#                                                                                                   

proc getA {Amin Amax} {

    # write current min and max values
    set fileID [open data1 w]
    puts $fileID "$Amin $Amax"
    close $fileID

    # invoke matlab
    if {[catch {exec matlab -nosplash -nodesktop -r "getA; quit"}]} {
	puts "Ignore this $msg"
    }

    # read matlab result
    set fileIN [open data2 r]
    set fileData [read $fileIN]
    close $fileIN

    #cleanup
    file delete data1
    file delete data2
    
    # return result
    set lines [split $fileData "\n"]
    set A [lindex $lines 0]  
}

set A [getA $Amin $Amax]

#                                                                                                   
# create the model & perform one analysis step                                                      
#                                                                                                   

model BasicBuilder -ndm 2 -ndf 2

node 1   0.0  0.0
node 2 144.0  0.0
node 3 168.0  0.0
node 4  72.0 96.0

fix 1 1 1
fix 2 1 1
fix 3 1 1
uniaxialMaterial Elastic 1 3000


element truss 1 1 4 $A 1
element truss 2 2 4 $A 1
element truss 3 3 4 $A 1

timeSeries Constant 1
pattern Plain 1 1 {
    load 4 100 -50
}

#create analysis

system BandSPD
numberer RCM
constraints Plain
integrator LoadControl 1.0
algorithm Linear
analysis Static

analyze 1

#                                                                                                   
# iterate until convergence is reached, or failure                                                  
#                                                                                                   

set xDisp [nodeDisp 4 1]
set done 0

while {$done == 0} {

    if {$xDisp > $targetDisp} {
	set Amin $A
    } else {
	set Amax $A
    }

    set diffA [expr $Amax-$Amin]
    if {$diffA  < $tol} {
	set done 2
      break
    }

    set A [getA $Amin $Amax]
    setParameter -value $A -ele A
    
    analyze 1
    set xDisp [nodeDisp 4 1]

    set diff [expr abs($xDisp - $targetDisp)]
    if {$diff < $tol} {
	set done 1
    }

    puts "$A $xDisp"
}

#                                                                                                   
# print success or failure                                                                          
#                  

if {$done == 1} {
    puts "SUCCESS: A = $A gives Displacement of $xDisp versus target of $targetDisp"
} else {
    puts "FAILURE: Hit limit A = $A gives Displacement of $xDisp versus target of $targetDisp"
}


The matlab scipt contained in the file getA.m is a simple divide-and-conquer.

load data1;
a=(data1(1,1)+data1(1,2))/2.0;
dlmwrite('data2', a);