OpenMDAO iteration loops and distributed architectures (ATC) in OpenMDAO

 I am attempting to implement ATC architecture into OpenMDAO v0.13.0. I've posted before about bugs and snags I've run into and have received help, but I'm still having issues. I cannot find any way to control iteration loops and the order in which sub problems are solved. I am trying to implement this MATLAB code into equivalent OpenMDAO: ``````% begin ATC iterations while conv == 0 % solve P0 [T11,T12,x0] = P0(x_old,R11,R12,v,w); % solve P11 [R11,y11,x11] = P11(x_old,T11,y12,v,w); % solve P12 [R12,y12,x12] = P12(x_old,T12,y11,v,w); % check convergence % new vector of optimization variables x_new = [T11 T12 x0 R11 y11 x11 R12 y12 x12]; % consistency constraints c = [T11-R11 T12-R12 y12-y11] c_norm(k+1) = norm(c); % choose criterion if c_norm(k+1) < ATC_eps % if (norm(x_new - x_old))/norm(x_old) < ATC_eps conv = 1; disp('convergence achieved') end if k > 1000 conv = 1; disp('maximum number of ATC iterations reached') end %prepare for reiteration x_old = x_new; %update penalty weights v = v + 2*w.*w.*c; w = beta*w; k = k + 1; end `````` Where P0(), P11(), and P12() are functions applying fconmin to the sub problems with the associated constraints. My code seems to be doing this process in reverse order, and is giving incorrect results. Any help you could provide on implementing distributed architectures, architectures with penalty weights, and decomposed problems in OpenMDAO would be greatly appreciated. In my code I used a top level assembly with a driver who's objective is to minimize the norm of c (the consistency constraints between targets and responses). Assembly File ``````from openmdao.lib.datatypes.api import Float, Array from openmdao.main.api import Assembly, set_as_top from openmdao.lib.drivers.api import SLSQPdriver import DisciplineTemplate from numpy import ones, array, float class Optimization0(Assembly): """my optimization template""" def configure(self): # Create Optimizer instance self.add('driver', SLSQPdriver()) self.driver.accuracy = 1.0e-6 self.add('temp', Float(1.0, iotype='in')) self.driver.add_parameter('temp', low=0., high=10.) # Create component instances self.add('Multipliers', DisciplineTemplate.Multipliers()) # Iteration Hierarchy self.driver.workflow.add(['Multipliers']) # SLSQP Flags self.driver.iprint = 0 class Optimization1(Assembly): """my optimization template""" v = Array(array([1.0,1.0,1.0]), dtype=float, desc = "Penalty", iotype='in') w = Array(array([1.0,1.0,1.0]), dtype=float, desc = "Penalty", iotype='in') pR11 = Float(1.0, desc='Target', iotype='in') pR12 = Float(1.0, desc='Coupling', iotype='in') def configure(self): # Create Optimizer instance self.add('driver', SLSQPdriver()) # Create component instances self.add('R01', DisciplineTemplate.R01()) self.add('R02', DisciplineTemplate.R02()) # Iteration Hierarchy self.driver.workflow.add(['R01', 'R02']) #x0[z4 z5 z7] self.driver.add_parameter('R01.z4', low=0., high=10.) self.driver.add_parameter(('R01.z5','R02.z5'), low = 0.0, high=10.0) self.driver.add_parameter('R02.z7', low=0., high=10.) #T11=z3 self.driver.add_parameter('R01.z3', low=0., high=10.) #T12=z6 self.driver.add_parameter('R02.z6', low=0., high=10.) #Targets #self.create_passthrough('R01.z3', 'pT11') #self.create_passthrough('R02.z6', 'pT12') self.add('pT11', Float(1.0, iotype='in')) self.connect('R01.z3', 'pT11') self.add('pT12', Float(1.0, iotype='in')) self.connect('R02.z6', 'pT12') #Constraints self.driver.add_constraint('(R01.z3)**(-2) + (R01.z4)**(2) - (R01.z5)**(2) <= 0') self.driver.add_constraint('(R02.z5)**(2) + (R02.z6)**(-2) - (R02.z7)**(2) <= 0') # SLSQP Flags self.driver.iprint = 0 class Optimization2(Assembly): """my optimization template""" v = Array(array([1.0,1.0,1.0]), dtype=float, desc = "Penalty", iotype='in') w = Array(array([1.0,1.0,1.0]), dtype=float, desc = "Penalty", iotype='in') pT11 = Float(1.0, desc='Target', iotype='in') py12 = Float(1.0, desc='Coupling', iotype='in') def configure(self): # Create Optimizer instance self.add('driver', SLSQPdriver()) # Create component instances self.add('R11', DisciplineTemplate.R11()) # Iteration Hierarchy self.driver.workflow.add(['R11']) #x11[z8,z9,z10] self.driver.add_parameter('R11.z8', low=0., high=10.) self.driver.add_parameter('R11.z9', low=0., high=10.) self.driver.add_parameter('R11.z10', low=0., high=10.) #y11=z11 self.driver.add_parameter('R11.z11', low=0., high=10.) #Out Targets #self.create_passthrough('R11.z3', 'pR11') self.add('pR11', Float(1.0, iotype='out')) self.connect('R11.z3', 'pR11') #self.create_passthrough('R11.z11', 'py11') self.add('py11', Float(1.0, iotype='out')) self.connect('R11.z11', 'py11') #Constraints self.driver.add_constraint('(R11.z8)**(2) + (R11.z9)**(2) - (R11.z11)**(2) <= 0') self.driver.add_constraint('(R11.z8)**(-2) + (R11.z10)**(2) - (R11.z11)**(2) <= 0') # SLSQP Flags self.driver.iprint = 0 class Optimization3(Assembly): """my optimization template""" v = Array(array([1.0,1.0,1.0]), dtype=float, desc = "Penalty", iotype='in') w = Array(array([1.0,1.0,1.0]), dtype=float, desc = "Penalty", iotype='in') pT12 = Float(1.0, desc='Target', iotype='in') py11 = Float(1.0, desc='Coupling', iotype='in') def configure(self): # Create Optimizer instance self.add('driver', SLSQPdriver()) # Create component instances self.add('R12', DisciplineTemplate.R12()) # Iteration Hierarchy self.driver.workflow.add(['R12']) #x12[z12,z13,z14] self.driver.add_parameter('R12.z12', low=0., high=10.) self.driver.add_parameter('R12.z13', low=0., high=10.) self.driver.add_parameter('R12.z14', low=0., high=10.) #y12=z11 self.driver.add_parameter('R12.z11', low=0., high=10.) #Targets #self.create_passthrough('R12.z6', 'pR12') self.add('pR12', Float(1.0, iotype='out')) self.connect('R12.z6', 'pR12') #self.create_passthrough('R12.z11', 'py12') self.add('py12', Float(1.0, iotype='out')) self.connect('R12.z11', 'py12') #Constraints self.driver.add_constraint('(R12.z11)**(2) + (R12.z12)**(-2) - (R12.z13)**(2) <= 0') self.driver.add_constraint('(R12.z11)**(2) + (R12.z12)**(2) - (R12.z14)**(2) <= 0') # SLSQP Flags self.driver.iprint = 0 if __name__ == "__main__": #Hierarchy Loop = set_as_top(Optimization0()) Loop.add('P0', Optimization1()) Loop.driver.workflow.add('P0') Loop.add('P11', Optimization2()) Loop.driver.workflow.add('P11') Loop.add('P12', Optimization3()) Loop.driver.workflow.add('P12') #connections Loop.connect('[P11.pT11-P11.pR11,P12.pT12-P12.pR12,P12.py12-P11.py11]','Multipliers.c') Loop.connect('Multipliers.v','P0.v') Loop.connect('Multipliers.w','P0.w') Loop.connect('P11.pR11','P0.pR11') Loop.connect('P12.pR12','P0.pR12') Loop.connect('Multipliers.v','P11.v') Loop.connect('Multipliers.w','P11.w') Loop.connect('P0.pT11','P11.pT11') Loop.connect('P12.py12','P11.py12') Loop.connect('Multipliers.v','P12.v') Loop.connect('Multipliers.w','P12.w') Loop.connect('P0.pT12','P12.pT12') Loop.connect('P11.py11','P12.py11') #Objectives Loop.driver.add_objective('(P0.pR11-P0.pT11)**2 + (P0.pR12-P0.pT12)**2 + (P12.py12-P11.py11)**2') Loop.P0.driver.add_objective('(R01.z1-0.0)**2 + (R02.z2-0.0)**2 + v[0]*(R01.z3-pR11) + v[1]*(R02.z6-pR12) + w[0]*w[0]*(R01.z3-pR11)**2 + w[1]*w[1]*(R02.z6-pR12)**2') Loop.P11.driver.add_objective('v[0]*(pT11-R11.z3) + v[2]*(py12-R11.z11) + w[0]*w[0]*(pT11-R11.z3)**2 + w[2]*w[2]*(py12-R11.z11)**2') Loop.P12.driver.add_objective('v[1]*(pT12-R12.z6) + v[2]*(R12.z11-py11) + w[1]*w[1]*(pT12-R12.z6)**2 + w[2]*w[2]*(R12.z11-py11)**2') import time tt = time.time() Loop.run() official=[2.80,3.03,2.35,0.76,0.87,2.79,0.95,0.97,0.87,0.80,1.30,0.84,1.75,1.54] z=ones(14) z[0]=Loop.P0.R01.z1 z[1]=Loop.P0.R02.z2 z[2]=Loop.P0.R01.z3 z[3]=Loop.P0.R01.z4 z[4]=Loop.P0.R01.z5 z[5]=Loop.P0.R02.z6 z[6]=Loop.P0.R02.z7 z[7]=Loop.P11.R11.z8 z[8]=Loop.P11.R11.z9 z[9]=Loop.P11.R11.z10 z[10]=Loop.P11.R11.z11 z[11]=Loop.P12.R12.z12 z[12]=Loop.P12.R12.z13 z[13]=Loop.P12.R12.z14 for i in range(0, 14): print "\t\t z%d & %f & %f \\\\" % (i+1, z[i],official[i]) print "\n" for i in range(0, 14): print "z%d=%f vs official value: %f" % (i+1, z[i], official[i]) print "\nElapsed time: ", time.time()-tt, "seconds \n" print "Function calls R01: %d, R02: %d, R11: %d, R12: %d"%(Loop.P0.R01.exec_count, Loop.P0.R02.exec_count, Loop.P11.R11.exec_count, Loop.P12.R12.exec_count) `````` Components File ``````from openmdao.main.api import Component from openmdao.lib.datatypes.api import Float, Array from numpy import ones, array, float class R01(Component): z3 = Float(1.0, iotype='in', desc='Coupling Design Variable') z4 = Float(1.0, iotype='in', desc='Local Design Variable') z5 = Float(1.0, iotype='in', desc='Shared Design Variable') z1 = Float(iotype='out', desc='Output of this Discipline') def execute(self): z3 = self.z3 z4 = self.z4 z5 = self.z5 self.z1 = ( z3**2.0 + z4**(-2.0) + z5**2.0 )**(1.0/2.0) print "R01" print "z3=", z3 print "z4=", z4 print "z5=", z5 print "z1=", self.z1 class R02(Component): z5 = Float(1.0, iotype='in', desc='Shared Design Variable') z6 = Float(1.0, iotype='in', desc='Coupling Design Variable') z7 = Float(1.0, iotype='in', desc='Local Design Variable') z2 = Float(iotype='out', desc='Output of this Discipline') def execute(self): z5 = self.z5 z6 = self.z6 z7 = self.z7 self.z2 = ( z5**2.0 + z6**(2.0) + z7**2.0 )**(1.0/2.0) print "R02" print "z5=", z5 print "z6=", z6 print "z7=", z7 print "z2=", self.z2 class R11(Component): z8 = Float(1.0, iotype='in', desc='Local Design Variable') z9 = Float(1.0, iotype='in', desc='Local Design Variable') z10 = Float(1.0, iotype='in', desc='Local Design Variable') z11 = Float(1.0, iotype='in', desc='Shared Design Variable') z3 = Float(iotype='out', desc='Output of this Discipline') def execute(self): z8 = self.z8 z9 = self.z9 z10 = self.z10 z11 = self.z11 self.z3 = ( z8**2.0 + z9**(-2.0) + z10**(-2.0) + z11**2.0 )**(1.0/2.0) print "R11" print "z8=", z8 print "z9=", z9 print "z10=", z10 print "z11=", z11 print "z3=", self.z3 class R12(Component): z11 = Float(1.0, iotype='in', desc='Shared Design Variable') z12 = Float(1.0, iotype='in', desc='Local Design Variable') z13 = Float(1.0, iotype='in', desc='Local Design Variable') z14 = Float(1.0, iotype='in', desc='Local Design Variable') z6 = Float(iotype='out', desc='Output of this Discipline') def execute(self): z11 = self.z11 z12 = self.z12 z13 = self.z13 z14 = self.z14 self.z6 = ( z11**2.0 + z12**2.0 + z13**2.0 + z14**2.0 )**(1.0/2.0) print "R12" print "z11=", z11 print "z12=", z12 print "z13=", z13 print "z14=", z14 print "z6=", self.z6 class Multipliers(Component): c = Array(array([1.0,1.0,1.0]), dtype=float, desc = "Residuals", iotype='in') v = Array(array([1.0,1.0,1.0]), dtype=float, desc = "Penalty", iotype='out') w = Array(array([1.0,1.0,1.0]), dtype=float, desc = "Penalty", iotype='out') beta = Float(1.25, desc='beta') #c = [T11-R11 T12-R12 y12-y11] def execute(self): v = self.v w = self.w c = self.c beta = self.beta self.v = v + 2*w*w*c self.w = beta*w print "Multipliers" print "c=", self.c print "v=", self.v print "w=", self.w `````` asked 06 Aug '15, 17:24 Grasshopper 37●2●9

 your ordering problem is related to your data passing. For instance you have `````` Loop.connect('P0.pT11','P11.pT11') `````` and Loop.connect('P11.pR11','P0.pR11') So you have a cycle between P0 and P1. OpenMDAO is probably ordering these in a different way then you would like. I think you want a jacobi like convergence behavior here, where each iteration just uses the values from the previous case and you loop around till it converges. You will need to look into using a the FixedPointIterator here. You can break one of your connections (the ones coming into P0) and then add the parameters and constraints to the FixedPointIterator accordingly. This should give you proper run order and actual convergence behavior. You can look at our sellar_bliss example to see a good way to set this up. answered 06 Aug '15, 20:00 justingray ♦♦ 1.8k●1●14

By Email:

Markdown Basics

• *italic* or _italic_
• **bold** or __bold__
• image?![alt text](/path/img.jpg "title")
• numbered list: 1. Foo 2. Bar
• to add a line break simply add two spaces to where you would like the new line to be.
• basic HTML tags are also supported

Tags:

Seen: 1,523 times

Last updated: 06 Aug '15, 20:00