source: issm/oecreview/Archive/13393-13976/ISSM-13899-13900.diff@ 14312

../trunk-jpl/test/NightlyRun/test3005.py

@@ -1 @@
-from MatlabFuncs import *
 from model import *
-from EnumDefinitions import *
-import numpy
 from triangle import *
 from setmask import *
 from parameterize import *
 from setflowequation import *
+from EnumDefinitions import *
 from solve import *
+from MatlabFuncs import *
 
-md=triangle(model(),'../Exp/Square.exp',150000)
+md=triangle(model(),'../Exp/Square.exp',150000.)
 md=setmask(md,'all','')
 md=parameterize(md,'../Par/SquareShelfConstrained.py')
 md=setflowequation(md,'macayeal','all')
 md.cluster=generic('name',oshostname(),'np',3)
-md.autodiff.isautodiff=true
+md.autodiff.isautodiff=True
 md=solve(md,PrognosticSolutionEnum())
 
-
-# Fields and tolerances to track changes
-
+#Fields and tolerances to track changes
 field_names     =['Thickness']
 field_tolerances=[1e-13]
 field_values=[\

../trunk-jpl/test/NightlyRun/test3019.py

@@ +1 @@
+import numpy
+from model import *
+from triangle import *
+from setmask import *
+from parameterize import *
+from setflowequation import *
+from independent import *
+from dependent import *
+from EnumDefinitions import *
+from solve import *
+from MatlabFuncs import *
+
+#test reverse scalar vs forward vectorial drivers in ADOLC, using the test3009 setup, equivalent to test109 setup.
+md=triangle(model(),'../Exp/Square.exp',100000.)
+md=setmask(md,'all','')
+md=parameterize(md,'../Par/SquareShelfConstrained.py')
+md=setflowequation(md,'macayeal','all')
+md.cluster=generic('name',oshostname(),'np',3)
+
+md.autodiff.isautodiff=True
+md.verbose.autodiff=True
+
+#first run scalar reverse mode: 
+md.autodiff.independents=[independent('name','Thickness','type','vertex','nods',md.mesh.numberofvertices)]
+md.autodiff.dependents=[dependent('name','MaxVel','type','scalar','fos_reverse_index',1)]
+md.autodiff.driver='fos_reverse'
+
+md=solve(md,TransientSolutionEnum())
+
+#recover jacobian: 
+jac_reverse=md.results['TransientSolution'][1]['AutodiffJacobian']
+
+#now run vectorial forward mode
+md.autodiff.independents=[independent('name','Thickness','type','vertex','nods',md.mesh.numberofvertices,'fov_forward_indices',numpy.arange(0,md.mesh.numberofvertices))]
+md.autodiff.dependents=[dependent('name','MaxVel','type','scalar')]
+md.autodiff.driver='fov_forward'
+
+md=solve(md,TransientSolutionEnum())
+
+#recover jacobian: 
+jac_forward=md.results['TransientSolution'][1]['AutodiffJacobian']
+
+#Fields and tolerances to track changes
+field_names     =['Jac Forward','Jac Reverse','Jac Forward - Reverse']
+field_tolerances=[1e-13,1e-13,1e-13]
+field_values=[jac_forward,jac_reverse,jac_forward-jac_reverse]

../trunk-jpl/test/NightlyRun/test3019.m

@@ -18 +18 @@
 %recover jacobian: 
 jac_reverse=md.results.TransientSolution(1).AutodiffJacobian;
 
-
 %now run vectorial forward mode
 md.autodiff.independents={independent('name','Thickness','type','vertex','nods',md.mesh.numberofvertices,'fov_forward_indices',(1:md.mesh.numberofvertices)')};
 md.autodiff.dependents={dependent('name','MaxVel','type','scalar')};
@@ -29 +28 @@
 %recover jacobian: 
 jac_forward=md.results.TransientSolution(1).AutodiffJacobian;
 
-
 %Fields and tolerances to track changes
 field_names     ={'Jac Forward','Jac Reverse','Jac Forward - Reverse'};
 field_tolerances={1e-13,1e-13,1e-13};

../trunk-jpl/test/NightlyRun/test3020.m

@@ -7 +7 @@
 md=parameterize(md,'../Par/SquareShelfConstrained.par');
 md=setflowequation(md,'macayeal','all');
 md.cluster=generic('name',oshostname(),'np',3);
-md.transient.requested_outputs=[IceVolumeEnum();MaxVelEnum];
+md.transient.requested_outputs=[IceVolumeEnum();MaxVelEnum()];
 md.verbose=verbose('autodiff',true);
-md.diagnostic.restol=.000001;
+md.diagnostic.restol=0.000001;
 
 %setup autodiff parameters
 index=1; %this is the scalar component we are checking against
@@ -23 +23 @@
         };
 md.autodiff.driver='fos_forward';
 
-%parameters for the step-wise devivative
-delta=.00001;
+%parameters for the step-wise derivative
+delta=0.00001;
 h1=md.geometry.thickness(index);
-h0=h1*(1-delta);
-h2=h1*(1+delta);
+h0=h1*(1.-delta);
+h2=h1*(1.+delta);
 deltaH=(h2-h0);
 
 %save model:
@@ -62 +62 @@
 dVdh_an=(V2-V0)/deltaH;
 dMaxVdh_an=(MaxV2-MaxV0)/deltaH;
 
-
 %evaluate derivative using ADOLC 
 md=md2;
 md.autodiff.isautodiff=true;
@@ -79 +78 @@
 disp(sprintf('dV/dh: analytical:  %16.16g\n       using adolc:  %16.16g\n',dVdh_an,dVdh_ad));
 disp(sprintf('dMaxV/dh: analytical:  %16.16g\n       using adolc:  %16.16g\n',dMaxVdh_an,dMaxVdh_ad));
 
-
 %Fields and tolerances to track changes
 field_names     ={'dV/dh-dV/dh0','dMaxV/dh-dMaxV/dh0'};
 field_tolerances={1e-13,1e-13};

../trunk-jpl/test/NightlyRun/test3006.py

@@ -1 @@
-from MatlabFuncs import *
 from model import *
-from EnumDefinitions import *
-import numpy
 from triangle import *
 from meshconvert import *
 from setmask import *
 from parameterize import *
 from setflowequation import *
+from EnumDefinitions import *
 from solve import *
+from MatlabFuncs import *
 
-md=triangle(model(),'../Exp/Square.exp',150000)
+md=triangle(model(),'../Exp/Square.exp',150000.)
 md=meshconvert(md)
 md=setmask(md,'all','')
 md=parameterize(md,'../Par/SquareShelfConstrained.py')
@@ -17 +16 @@
 md.cluster=generic('name',oshostname(),'np',3)
 md.prognostic.stabilization=3
 md.prognostic.spcthickness=md.geometry.thickness
-md.autodiff.isautodiff=true
+md.autodiff.isautodiff=True
 md=solve(md,PrognosticSolutionEnum())
 
-
-# Fields and tolerances to track changes
-
+#Fields and tolerances to track changes
 field_names     =['Thickness']
 field_tolerances=[1e-13]
 field_values=[\

../trunk-jpl/test/NightlyRun/test3001.py

@@ -1 @@
-from MatlabFuncs import *
 from model import *
-from EnumDefinitions import *
-import numpy
 from triangle import *
 from setmask import *
 from parameterize import *
 from setflowequation import *
+from EnumDefinitions import *
 from solve import *
+from MatlabFuncs import *
 
-md=triangle(model(),'../Exp/Square.exp',50000)
+md=triangle(model(),'../Exp/Square.exp',50000.)
 md=setmask(md,'all','')
 md=parameterize(md,'../Par/SquareShelfConstrained.py')
 md=setflowequation(md,'macayeal','all')
 md.cluster=generic('name',oshostname(),'np',3)
 md.diagnostic.requested_outputs=StressTensorEnum()
-md.autodiff.isautodiff=true
+md.autodiff.isautodiff=True
 md=solve(md,DiagnosticSolutionEnum())
 
-
-# Fields and tolerances to track changes
-
+#Fields and tolerances to track changes
 field_names     =['Vx','Vy','Vel','Pressure',\
         'StressTensorxx','StressTensoryy','StressTensorxy']
 field_tolerances=[1e-13,1e-13,1e-13,1e-13,\

../trunk-jpl/test/NightlyRun/test3015.py

@@ +1 @@
+import numpy
+import copy
+from model import *
+from triangle import *
+from setmask import *
+from parameterize import *
+from setflowequation import *
+from independent import *
+from dependent import *
+from SetIceShelfBC import *
+from EnumDefinitions import *
+from solve import *
+from MatlabFuncs import *
+
+#This test runs test3005 with autodiff on, and checks that 
+#the value of the scalar forward difference match a step-wise differential
+
+#First configure
+md=triangle(model(),'../Exp/Square.exp',50000.)
+md=setmask(md,'all','')
+md=parameterize(md,'../Par/SquareShelfConstrained.py')
+md=setflowequation(md,'macayeal','all')
+md.cluster=generic('name',oshostname(),'np',3)
+md.prognostic.requested_outputs=IceVolumeEnum()
+md.verbose=verbose('autodiff',True)
+
+#setup autodiff parameters
+index=0 #this is the scalar component we are checking against
+md.autodiff.independents=[\
+        independent('name','Thickness','type','vertex','nods',md.mesh.numberofvertices,'fos_forward_index',index)
+        ]
+
+md.autodiff.dependents=[\
+        dependent('name','IceVolume','type','scalar')\
+        ]
+md.autodiff.driver='fos_forward'
+
+#parameters for the step-wise derivative
+delta=0.001
+h1=md.geometry.thickness[index]
+h0=h1*(1.-delta)
+h2=h1*(1.+delta)
+deltaH=(h2-h0)
+
+#save model
+md2=copy.deepcopy(md)
+
+#evaluate derivative by forward and backward stepping 
+#forward
+md=copy.deepcopy(md2)
+md.autodiff.isautodiff=False
+md.geometry.thickness[index]=h0
+md.geometry.bed=-md.materials.rho_ice/md.materials.rho_water*md.geometry.thickness
+md.geometry.surface=md.geometry.bed+md.geometry.thickness
+md=SetIceShelfBC(md)
+
+md=solve(md,PrognosticSolutionEnum())
+V0=md.results['PrognosticSolution'][1]['IceVolume']
+
+#backward
+md=copy.deepcopy(md2)
+md.autodiff.isautodiff=False
+md.geometry.thickness[index]=h2
+md.geometry.bed=-md.materials.rho_ice/md.materials.rho_water*md.geometry.thickness
+md.geometry.surface=md.geometry.bed+md.geometry.thickness
+md=SetIceShelfBC(md)
+
+md=solve(md,PrognosticSolutionEnum())
+V2=md.results['PrognosticSolution'][1]['IceVolume']
+
+#compute resulting derivative
+dVdh_an=(V2-V0)/deltaH
+
+#evaluate derivative using ADOLC 
+md=md2
+md.autodiff.isautodiff=True
+md.geometry.thickness[index]=h1
+md.geometry.bed=-md.materials.rho_ice/md.materials.rho_water*md.geometry.thickness
+md.geometry.surface=md.geometry.bed+md.geometry.thickness
+md=SetIceShelfBC(md)
+
+md=solve(md,PrognosticSolutionEnum())
+#retrieve directly
+dVdh_ad=md.results['PrognosticSolution'][1]['AutodiffJacobian']
+
+print "dV/dh: analytical:  #16.16g\n       using adolc:  #16.16g\n" % (dVdh_an,dVdh_ad)
+
+#Fields and tolerances to track changes
+field_names     =['dV/dh-dV/dh0']
+field_tolerances=[1e-13]
+field_values=[dVdh_ad-dVdh_an]

../trunk-jpl/test/NightlyRun/test3010.py

@@ -1 @@
-from MatlabFuncs import *
 from model import *
-from EnumDefinitions import *
-import numpy
 from triangle import *
 from setmask import *
 from parameterize import *
 from setflowequation import *
+from EnumDefinitions import *
 from solve import *
+from MatlabFuncs import *
 
-md=triangle(model(),'../Exp/Square.exp',150000)
+md=triangle(model(),'../Exp/Square.exp',150000.)
 md=setmask(md,'all','')
 md=parameterize(md,'../Par/SquareShelfConstrained.py')
 md=setflowequation(md,'macayeal','all')
 md.cluster=generic('name',oshostname(),'np',3)
 md.transient.requested_outputs=IceVolumeEnum()
 
-
-md.autodiff.isautodiff=true
+md.autodiff.isautodiff=True
 md=solve(md,TransientSolutionEnum())
 
-
-# Fields and tolerances to track changes
-
+#Fields and tolerances to track changes
 field_names     =['Vx1','Vy1','Vel1','Pressure1','Bed1','Surface1','Thickness1','Volume1','Vx2','Vy2','Vel2','Pressure2','Bed2','Surface2','Thickness2','Volume2','Vx3','Vy3','Vel3','Pressure3','Bed3','Surface3','Thickness3','Volume3']
 field_tolerances=[1e-13,1e-13,1e-13,1e-13,1e-13,1e-13,1e-13,1e-13,\
                                                 1e-13,1e-13,1e-13,1e-13,1e-13,1e-13,1e-13,1e-13,\

../trunk-jpl/test/NightlyRun/test3007.py

@@ -1 @@
-from MatlabFuncs import *
 from model import *
-from EnumDefinitions import *
-import numpy
 from triangle import *
 from setmask import *
 from parameterize import *
 from setflowequation import *
+from EnumDefinitions import *
 from solve import *
+from MatlabFuncs import *
 
-md=triangle(model(),'../Exp/Square.exp',150000)
+md=triangle(model(),'../Exp/Square.exp',150000.)
 md=setmask(md,'all','')
 md=parameterize(md,'../Par/SquareShelfConstrained.py')
 md=setflowequation(md,'macayeal','all')
 md.extrude(5,3.)
 md.cluster=generic('name',oshostname(),'np',3)
-md.autodiff.isautodiff=true
+md.autodiff.isautodiff=True
 md=solve(md,PrognosticSolutionEnum())
 
-
-# Fields and tolerances to track changes
-
+#Fields and tolerances to track changes
 field_names     =['Thickness']
 field_tolerances=[1e-13]
 field_values=[\

../trunk-jpl/test/NightlyRun/test3002.py

@@ -1 @@
-from MatlabFuncs import *
 from model import *
-from EnumDefinitions import *
-import numpy
 from triangle import *
 from setmask import *
 from parameterize import *
 from setflowequation import *
+from EnumDefinitions import *
 from solve import *
+from MatlabFuncs import *
 
-md=triangle(model(),'../Exp/Square.exp',180000)
+md=triangle(model(),'../Exp/Square.exp',180000.)
 md=setmask(md,'all','')
 md=parameterize(md,'../Par/SquareShelfConstrained.py')
 md.extrude(3,2.)
 md=setflowequation(md,'macayeal','all')
 md.cluster=generic('name',oshostname(),'np',3)
-md.autodiff.isautodiff=true
+md.autodiff.isautodiff=True
 md=solve(md,DiagnosticSolutionEnum())
 
-
-# Fields and tolerances to track changes
-
+#Fields and tolerances to track changes
 field_names     =['Vx','Vy','Vz','Vel','Pressure']
 field_tolerances=[1e-13,1e-13,1e-13,1e-13,1e-13]
 field_values=[\

../trunk-jpl/test/NightlyRun/test3020.py

@@ +1 @@
+import numpy
+import copy
+from model import *
+from triangle import *
+from setmask import *
+from parameterize import *
+from setflowequation import *
+from independent import *
+from dependent import *
+from SetIceShelfBC import *
+from EnumDefinitions import *
+from solve import *
+from MatlabFuncs import *
+
+#This test runs test3020 with autodiff on, and checks that 
+#the value of the scalar forward difference match a step-wise differential
+
+#First configure
+md=triangle(model(),'../Exp/Square.exp',150000.)
+md=setmask(md,'all','')
+md=parameterize(md,'../Par/SquareShelfConstrained.py')
+md=setflowequation(md,'macayeal','all')
+md.cluster=generic('name',oshostname(),'np',3)
+md.transient.requested_outputs=[IceVolumeEnum(),MaxVelEnum()]
+md.verbose=verbose('autodiff',True)
+md.diagnostic.restol=0.000001
+
+#setup autodiff parameters
+index=0 #this is the scalar component we are checking against
+md.autodiff.independents=[\
+        independent('name','Thickness','type','vertex','nods',md.mesh.numberofvertices,'fos_forward_index',index)
+        ]
+
+md.autodiff.dependents=[\
+        dependent('name','IceVolume','type','scalar'),\
+        dependent('name','MaxVel','type','scalar')\
+        ]
+md.autodiff.driver='fos_forward'
+
+#parameters for the step-wise derivative
+delta=0.00001
+h1=md.geometry.thickness[index]
+h0=h1*(1.-delta)
+h2=h1*(1.+delta)
+deltaH=(h2-h0)
+
+#save model:
+md2=copy.deepcopy(md)
+
+#evaluate derivative by forward and backward stepping 
+#forward
+md=copy.deepcopy(md2)
+md.autodiff.isautodiff=False
+md.geometry.thickness[index]=h0
+md.geometry.bed=-md.materials.rho_ice/md.materials.rho_water*md.geometry.thickness
+md.geometry.surface=md.geometry.bed+md.geometry.thickness
+md=SetIceShelfBC(md)
+
+md=solve(md,TransientSolutionEnum())
+V0=md.results['TransientSolution'][3]['IceVolume']
+MaxV0=md.results['TransientSolution'][3]['MaxVel']
+
+#backward
+md=copy.deepcopy(md2)
+md.autodiff.isautodiff=False
+md.geometry.thickness[index]=h2
+md.geometry.bed=-md.materials.rho_ice/md.materials.rho_water*md.geometry.thickness
+md.geometry.surface=md.geometry.bed+md.geometry.thickness
+md=SetIceShelfBC(md)
+
+md=solve(md,TransientSolutionEnum())
+V2=md.results['TransientSolution'][3]['IceVolume']
+MaxV2=md.results['TransientSolution'][3]['MaxVel']
+
+#compute resulting derivative
+dVdh_an=(V2-V0)/deltaH
+dMaxVdh_an=(MaxV2-MaxV0)/deltaH
+
+#evaluate derivative using ADOLC 
+md=copy.deepcopy(md2)
+md.autodiff.isautodiff=True
+md.geometry.thickness[index]=h1
+md.geometry.bed=-md.materials.rho_ice/md.materials.rho_water*md.geometry.thickness
+md.geometry.surface=md.geometry.bed+md.geometry.thickness
+md=SetIceShelfBC(md)
+
+md=solve(md,TransientSolutionEnum())
+#retrieve directly
+dVdh_ad=md.results['TransientSolution'][1]['AutodiffJacobian'][0]
+dMaxVdh_ad=md.results['TransientSolution'][1]['AutodiffJacobian'][1]
+
+print "dV/dh: analytical:  %16.16g\n       using adolc:  %16.16g\n" % (dVdh_an,dVdh_ad)
+print "dMaxV/dh: analytical:  %16.16g\n       using adolc:  %16.16g\n" % (dMaxVdh_an,dMaxVdh_ad)
+
+#Fields and tolerances to track changes
+field_names     =['dV/dh-dV/dh0','dMaxV/dh-dMaxV/dh0']
+field_tolerances=[1e-13,1e-13]
+field_values=[dVdh_ad-dVdh_an,dMaxVdh_an-dMaxVdh_ad]

../trunk-jpl/test/NightlyRun/test3008.py

@@ -1 @@
-from MatlabFuncs import *
 from model import *
-from EnumDefinitions import *
-import numpy
 from triangle import *
 from setmask import *
 from parameterize import *
 from setflowequation import *
+from EnumDefinitions import *
 from solve import *
+from MatlabFuncs import *
 
-md=triangle(model(),'../Exp/Square.exp',180000)
+md=triangle(model(),'../Exp/Square.exp',180000.)
 md=setmask(md,'all','')
 md=parameterize(md,'../Par/SquareShelfConstrained.py')
 md.extrude(3,1.)
 md=setflowequation(md,'macayeal','all')
 md.timestepping.time_step=0
 md.cluster=generic('name',oshostname(),'np',3)
-md.autodiff.isautodiff=true
+md.autodiff.isautodiff=True
 md=solve(md,ThermalSolutionEnum())
 
-
-# Fields and tolerances to track changes
-
+#Fields and tolerances to track changes
 field_names     =['Temperature','BasalforcingsMeltingRate']
 field_tolerances=[1e-13,1e-13]
 field_values=[\

../trunk-jpl/test/NightlyRun/test3015.m

@@ -21 +21 @@
         };
 md.autodiff.driver='fos_forward';
 
-%parameters for the step-wise devivative
-delta=.001;
+%parameters for the step-wise derivative
+delta=0.001;
 h1=md.geometry.thickness(index);
-h0=h1*(1-delta);
-h2=h1*(1+delta);
+h0=h1*(1.-delta);
+h2=h1*(1.+delta);
 deltaH=(h2-h0);
 
 %save model:

../trunk-jpl/test/NightlyRun/test3003.py

@@ -1 @@
-from MatlabFuncs import *
 from model import *
-from EnumDefinitions import *
-import numpy
 from triangle import *
 from setmask import *
 from parameterize import *
 from setflowequation import *
+from EnumDefinitions import *
 from solve import *
+from MatlabFuncs import *
 
-md=triangle(model(),'../Exp/Square.exp',180000)
+md=triangle(model(),'../Exp/Square.exp',180000.)
 md=setmask(md,'all','')
 md=parameterize(md,'../Par/SquareShelfConstrained.py')
 md.extrude(3,2.)
 md=setflowequation(md,'pattyn','all')
 md.cluster=generic('name',oshostname(),'np',3)
 md.diagnostic.requested_outputs=StressTensorEnum()
-md.autodiff.isautodiff=true
+md.autodiff.isautodiff=True
 md=solve(md,DiagnosticSolutionEnum())
 
-
-# Fields and tolerances to track changes
-
+#Fields and tolerances to track changes
 field_names     =['Vx','Vy','Vz','Vel','Pressure',\
         'StressTensorxx','StressTensoryy','StressTensorzz','StressTensorxy','StressTensorxz','StressTensoryz']
 field_tolerances=[1e-09,1e-09,1e-09,1e-09,1e-09,\

../trunk-jpl/test/NightlyRun/test3009.py

@@ -1 @@
-from MatlabFuncs import *
 from model import *
-from EnumDefinitions import *
-import numpy
 from triangle import *
 from setmask import *
 from parameterize import *
 from setflowequation import *
+from EnumDefinitions import *
 from solve import *
+from MatlabFuncs import *
 
-md=triangle(model(),'../Exp/Square.exp',180000)
+md=triangle(model(),'../Exp/Square.exp',180000.)
 md=setmask(md,'all','')
 md=parameterize(md,'../Par/SquareShelfConstrained.py')
 md.extrude(3,1.)
@@ -18 +17 @@
 md.transient.isprognostic=0
 md.transient.isthermal=1
 md.transient.isgroundingline=0
-md.autodiff.isautodiff=true
+md.autodiff.isautodiff=True
 md=solve(md,TransientSolutionEnum())
 
-
-# Fields and tolerances to track changes
-
+#Fields and tolerances to track changes
 field_names     =['Temperature','BasalforcingsMeltingRate']
 field_tolerances=[1e-13,1e-13]
 field_values=[\

../trunk-jpl/test/NightlyRun/test3004.py

@@ -1 @@
-from MatlabFuncs import *
 from model import *
-from EnumDefinitions import *
-import numpy
 from triangle import *
 from setmask import *
 from parameterize import *
 from setflowequation import *
+from EnumDefinitions import *
 from solve import *
+from MatlabFuncs import *
 
-md=triangle(model(),'../Exp/Square.exp',180000)
+md=triangle(model(),'../Exp/Square.exp',180000.)
 md=setmask(md,'all','')
 md=parameterize(md,'../Par/SquareShelfConstrained.py')
 md.extrude(3,2.)
 md=setflowequation(md,'stokes','all')
 md.cluster=generic('name',oshostname(),'np',3)
-md.autodiff.isautodiff=true
+md.autodiff.isautodiff=True
 md=solve(md,DiagnosticSolutionEnum())
 
-
-# Fields and tolerances to track changes
-
+#Fields and tolerances to track changes
 field_names     =['Vx','Vy','Vz','Vel','Pressure']
 field_tolerances=[1e-08,1e-08,1e-07,1e-08,1e-08]
 field_values=[\