Index: /issm/trunk-jpl/src/c/analyses/AdjointBalancethickness2Analysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/AdjointBalancethickness2Analysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/AdjointBalancethickness2Analysis.cpp (revision 25317) @@ -116,5 +116,5 @@ _error_("not implemented yet"); }/*}}}*/ -void AdjointBalancethickness2Analysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void AdjointBalancethickness2Analysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ /*The gradient of the cost function is calculated in 2 parts. * Index: /issm/trunk-jpl/src/c/analyses/AdjointBalancethickness2Analysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/AdjointBalancethickness2Analysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/AdjointBalancethickness2Analysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void GradientJdHdt(Element* element,Vector* gradient,int control_index); void GradientJOmega(Element* element,Vector* gradient,int control_index); Index: /issm/trunk-jpl/src/c/analyses/AdjointBalancethicknessAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/AdjointBalancethicknessAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/AdjointBalancethicknessAnalysis.cpp (revision 25317) @@ -156,5 +156,5 @@ _error_("not implemented yet"); }/*}}}*/ -void AdjointBalancethicknessAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void AdjointBalancethicknessAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ /*The gradient of the cost function is calculated in 2 parts. * @@ -194,8 +194,8 @@ /*Deal with second term*/ switch(control_type){ - case BalancethicknessSpcthicknessEnum: GradientJDirichlet(element,gradient,control_index); break; - case BalancethicknessThickeningRateEnum: GradientJDhDt(element,gradient,control_index); break; - case VxEnum: GradientJVx( element,gradient,control_index); break; - case VyEnum: GradientJVy( element,gradient,control_index); break; + case BalancethicknessSpcthicknessEnum: GradientJDirichlet(element,gradient,control_interp,control_index); break; + case BalancethicknessThickeningRateEnum: GradientJDhDt(element,gradient,control_interp,control_index); break; + case VxEnum: GradientJVx( element,gradient,control_interp,control_index); break; + case VyEnum: GradientJVy( element,gradient,control_interp,control_index); break; default: _error_("control type not supported yet: " << EnumToStringx(control_type)); } @@ -205,5 +205,5 @@ }/*}}}*/ -void AdjointBalancethicknessAnalysis::GradientJDirichlet(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointBalancethicknessAnalysis::GradientJDirichlet(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Fetch number of vertices for this finite element*/ @@ -242,5 +242,5 @@ delete Ke; }/*}}}*/ -void AdjointBalancethicknessAnalysis::GradientJDhDt(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointBalancethicknessAnalysis::GradientJDhDt(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Fetch number of vertices for this finite element*/ @@ -266,5 +266,5 @@ xDelete(vertexpidlist); }/*}}}*/ -void AdjointBalancethicknessAnalysis::GradientJVx(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointBalancethicknessAnalysis::GradientJVx(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Intermediaries*/ @@ -313,5 +313,5 @@ delete gauss; }/*}}}*/ -void AdjointBalancethicknessAnalysis::GradientJVy(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointBalancethicknessAnalysis::GradientJVy(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Intermediaries*/ Index: /issm/trunk-jpl/src/c/analyses/AdjointBalancethicknessAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/AdjointBalancethicknessAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/AdjointBalancethicknessAnalysis.h (revision 25317) @@ -27,9 +27,9 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); - void GradientJDirichlet(Element* element,Vector* gradient,int control_index); - void GradientJDhDt(Element* element,Vector* gradient,int control_index); - void GradientJVx(Element* element,Vector* gradient,int control_index); - void GradientJVy(Element* element,Vector* gradient,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); + void GradientJDirichlet(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDhDt(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJVx(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJVy(Element* element,Vector* gradient,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/AdjointHorizAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/AdjointHorizAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/AdjointHorizAnalysis.cpp (revision 25317) @@ -1078,5 +1078,5 @@ _error_("not implemented yet"); }/*}}}*/ -void AdjointHorizAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ /*The gradient of the cost function is calculated in 2 parts. * @@ -1103,6 +1103,6 @@ if(control_type!=MaterialsRheologyBbarEnum && control_type!=FrictionCoefficientEnum && - control_type!=FrictionCEnum && - control_type!=FrictionAsEnum && + control_type!=FrictionCEnum && + control_type!=FrictionAsEnum && control_type!=DamageDbarEnum && control_type!=MaterialsRheologyBEnum){ @@ -1117,8 +1117,8 @@ case SurfaceLogVxVyMisfitEnum: /*Nothing, \partial J/\partial k = 0*/ break; case SurfaceAverageVelMisfitEnum: /*Nothing, \partial J/\partial k = 0*/ break; - case DragCoefficientAbsGradientEnum: GradientJDragGradient(element,gradient,control_index); break; - case RheologyBbarAbsGradientEnum: GradientJBbarGradient(element,gradient,control_index); break; - case RheologyBAbsGradientEnum: GradientJBGradient(element,gradient,control_index); break; - case RheologyBInitialguessMisfitEnum: GradientJBinitial(element,gradient,control_index); break; + case DragCoefficientAbsGradientEnum: GradientJDragGradient(element,gradient,control_interp,control_index); break; + case RheologyBbarAbsGradientEnum: GradientJBbarGradient(element,gradient,control_interp,control_index); break; + case RheologyBAbsGradientEnum: GradientJBGradient(element,gradient,control_interp,control_index); break; + case RheologyBInitialguessMisfitEnum: GradientJBinitial(element,gradient,control_interp,control_index); break; default: _error_("response " << EnumToStringx(responses[resp]) << " not supported yet"); } @@ -1129,8 +1129,8 @@ case FrictionCEnum: switch(approximation){ - case SSAApproximationEnum: GradientJDragSSA(element,gradient,control_index); break; - case L1L2ApproximationEnum:GradientJDragL1L2(element,gradient,control_index); break; - case HOApproximationEnum: GradientJDragHO( element,gradient,control_index); break; - case FSApproximationEnum: GradientJDragFS( element,gradient,control_index); break; + case SSAApproximationEnum: GradientJDragSSA(element,gradient,control_interp,control_index); break; + case L1L2ApproximationEnum:GradientJDragL1L2(element,gradient,control_interp,control_index); break; + case HOApproximationEnum: GradientJDragHO( element,gradient,control_interp,control_index); break; + case FSApproximationEnum: GradientJDragFS( element,gradient,control_interp,control_index); break; case NoneApproximationEnum: /*Gradient is 0*/ break; default: _error_("approximation " << EnumToStringx(approximation) << " not supported yet"); @@ -1139,8 +1139,8 @@ case FrictionAsEnum: switch(approximation){ - case SSAApproximationEnum: GradientJDragHydroSSA(element,gradient,control_index); break; - case L1L2ApproximationEnum:GradientJDragHydroL1L2(element,gradient,control_index); break; - case HOApproximationEnum: GradientJDragHydroHO( element,gradient,control_index); break; - case FSApproximationEnum: GradientJDragHydroFS( element,gradient,control_index); break; + case SSAApproximationEnum: GradientJDragHydroSSA(element,gradient,control_interp,control_index); break; + case L1L2ApproximationEnum:GradientJDragHydroL1L2(element,gradient,control_interp,control_index); break; + case HOApproximationEnum: GradientJDragHydroHO( element,gradient,control_interp,control_index); break; + case FSApproximationEnum: GradientJDragHydroFS( element,gradient,control_interp,control_index); break; case NoneApproximationEnum: /*Gradient is 0*/ break; default: _error_("approximation " << EnumToStringx(approximation) << " not supported yet"); @@ -1149,8 +1149,8 @@ case MaterialsRheologyBbarEnum: switch(approximation){ - case SSAApproximationEnum: GradientJBbarSSA(element,gradient,control_index); break; - case L1L2ApproximationEnum:GradientJBbarL1L2(element,gradient,control_index); break; - case HOApproximationEnum: GradientJBbarHO( element,gradient,control_index); break; - case FSApproximationEnum: GradientJBbarFS( element,gradient,control_index); break; + case SSAApproximationEnum: GradientJBbarSSA(element,gradient,control_interp,control_index); break; + case L1L2ApproximationEnum:GradientJBbarL1L2(element,gradient,control_interp,control_index); break; + case HOApproximationEnum: GradientJBbarHO( element,gradient,control_interp,control_index); break; + case FSApproximationEnum: GradientJBbarFS( element,gradient,control_interp,control_index); break; case NoneApproximationEnum: /*Gradient is 0*/ break; default: _error_("approximation " << EnumToStringx(approximation) << " not supported yet"); @@ -1159,7 +1159,7 @@ case MaterialsRheologyBEnum: switch(approximation){ - case SSAApproximationEnum: GradientJBSSA(element,gradient,control_index); break; - case HOApproximationEnum: GradientJBHO( element,gradient,control_index); break; - case FSApproximationEnum: GradientJBFS( element,gradient,control_index); break; + case SSAApproximationEnum: GradientJBSSA(element,gradient,control_interp,control_index); break; + case HOApproximationEnum: GradientJBHO( element,gradient,control_interp,control_index); break; + case FSApproximationEnum: GradientJBFS( element,gradient,control_interp,control_index); break; case NoneApproximationEnum: /*Gradient is 0*/ break; default: _error_("approximation " << EnumToStringx(approximation) << " not supported yet"); @@ -1168,5 +1168,5 @@ case DamageDbarEnum: switch(approximation){ - case SSAApproximationEnum: GradientJDSSA(element,gradient,control_index); break; + case SSAApproximationEnum: GradientJDSSA(element,gradient,control_interp,control_index); break; case NoneApproximationEnum: /*Gradient is 0*/ break; default: _error_("approximation " << EnumToStringx(approximation) << " not supported yet"); @@ -1180,9 +1180,9 @@ }/*}}}*/ -void AdjointHorizAnalysis::GradientJBbarFS(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJBbarFS(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*WARNING: We use SSA as an estimate for now*/ - this->GradientJBbarSSA(element,gradient,control_index); -}/*}}}*/ -void AdjointHorizAnalysis::GradientJBbarGradient(Element* element,Vector* gradient,int control_index){/*{{{*/ + this->GradientJBbarSSA(element,gradient,control_interp,control_index); +}/*}}}*/ +void AdjointHorizAnalysis::GradientJBbarGradient(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Intermediaries*/ @@ -1259,15 +1259,15 @@ if(domaintype!=Domain2DhorizontalEnum){basalelement->DeleteMaterials(); delete basalelement;}; }/*}}}*/ -void AdjointHorizAnalysis::GradientJBbarL1L2(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJBbarL1L2(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Same as SSA*/ - return this->GradientJBbarSSA(element,gradient,control_index); -}/*}}}*/ -void AdjointHorizAnalysis::GradientJBbarHO(Element* element,Vector* gradient,int control_index){/*{{{*/ + return this->GradientJBbarSSA(element,gradient,control_interp,control_index); +}/*}}}*/ +void AdjointHorizAnalysis::GradientJBbarHO(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*WARNING: We use SSA as an estimate for now*/ - this->GradientJBbarSSA(element,gradient,control_index); -}/*}}}*/ -void AdjointHorizAnalysis::GradientJBbarSSA(Element* element,Vector* gradient,int control_index){/*{{{*/ + this->GradientJBbarSSA(element,gradient,control_interp,control_index); +}/*}}}*/ +void AdjointHorizAnalysis::GradientJBbarSSA(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Intermediaries*/ @@ -1336,12 +1336,31 @@ /*Build gradient vector (actually -dJ/dB): */ - for(int i=0;iweight*basis[i]; + _assert_(!xIsNan(ge[i])); + } + } + else if(control_interp==P0Enum){ + ge[0]+=-dmudB*thickness*( (2*dvx[0]+dvy[1])*2*dadjx[0]+(dvx[1]+dvy[0])*(dadjx[1]+dadjy[0])+(2*dvy[1]+dvx[0])*2*dadjy[1] - )*Jdet*gauss->weight*basis[i]; - _assert_(!xIsNan(ge[i])); - } - } - gradient->SetValues(numvertices,vertexpidlist,ge,ADD_VAL); + )*Jdet*gauss->weight; + _assert_(!xIsNan(ge[0])); + } + else{ + _error_("not supported"); + } + } + if(control_interp==P1Enum){ + gradient->SetValues(numvertices,vertexpidlist,ge,ADD_VAL); + } + else if(control_interp==P0Enum){ + gradient->SetValue(vertexpidlist[0],ge[0],ADD_VAL); + } + else{ + _error_("not supported"); + } /*Clean up and return*/ @@ -1353,9 +1372,9 @@ if(domaintype!=Domain2DhorizontalEnum){basalelement->DeleteMaterials(); delete basalelement;}; }/*}}}*/ -void AdjointHorizAnalysis::GradientJBFS(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJBFS(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*WARNING: We use HO as an estimate for now*/ - this->GradientJBHO(element,gradient,control_index); -}/*}}}*/ -void AdjointHorizAnalysis::GradientJBGradient(Element* element,Vector* gradient,int control_index){/*{{{*/ + this->GradientJBHO(element,gradient,control_interp,control_index); +}/*}}}*/ +void AdjointHorizAnalysis::GradientJBGradient(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Intermediaries*/ @@ -1423,5 +1442,5 @@ delete gauss; }/*}}}*/ -void AdjointHorizAnalysis::GradientJBinitial(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJBinitial(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Intermediaries*/ @@ -1487,5 +1506,5 @@ delete gauss; }/*}}}*/ -void AdjointHorizAnalysis::GradientJBHO(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJBHO(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Intermediaries*/ int domaintype,dim; @@ -1563,5 +1582,5 @@ delete gauss; }/*}}}*/ -void AdjointHorizAnalysis::GradientJBSSA(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJBSSA(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Intermediaries*/ @@ -1647,5 +1666,5 @@ if(domaintype!=Domain2DhorizontalEnum){basalelement->DeleteMaterials(); delete basalelement;}; }/*}}}*/ -void AdjointHorizAnalysis::GradientJDragGradient(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJDragGradient(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*return if floating (gradient is 0)*/ @@ -1741,5 +1760,5 @@ }/*}}}*/ -void AdjointHorizAnalysis::GradientJDragFS(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJDragFS(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*return if floating or not on bed (gradient is 0)*/ @@ -1837,10 +1856,10 @@ delete friction; }/*}}}*/ -void AdjointHorizAnalysis::GradientJDragL1L2(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJDragL1L2(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Same as SSA*/ - return this->GradientJDragSSA(element,gradient,control_index); -}/*}}}*/ -void AdjointHorizAnalysis::GradientJDragHO(Element* element,Vector* gradient,int control_index){/*{{{*/ + return this->GradientJDragSSA(element,gradient,control_interp,control_index); +}/*}}}*/ +void AdjointHorizAnalysis::GradientJDragHO(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*return if floating or not on bed (gradient is 0)*/ @@ -1916,5 +1935,5 @@ delete friction; }/*}}}*/ -void AdjointHorizAnalysis::GradientJDragSSA(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJDragSSA(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*return if floating (gradient is 0)*/ @@ -2018,6 +2037,5 @@ if(domaintype!=Domain2DhorizontalEnum){basalelement->DeleteMaterials(); delete basalelement;}; }/*}}}*/ - -void AdjointHorizAnalysis::GradientJDragHydroFS(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJDragHydroFS(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*return if floating or not on bed (gradient is 0)*/ @@ -2113,10 +2131,10 @@ delete friction; }/*}}}*/ -void AdjointHorizAnalysis::GradientJDragHydroL1L2(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJDragHydroL1L2(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Same as SSA*/ - return this->GradientJDragSSA(element,gradient,control_index); -}/*}}}*/ -void AdjointHorizAnalysis::GradientJDragHydroHO(Element* element,Vector* gradient,int control_index){/*{{{*/ + return this->GradientJDragSSA(element,gradient,control_interp,control_index); +}/*}}}*/ +void AdjointHorizAnalysis::GradientJDragHydroHO(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*return if floating or not on bed (gradient is 0)*/ @@ -2190,6 +2208,5 @@ delete friction; }/*}}}*/ - -void AdjointHorizAnalysis::GradientJDragHydroSSA(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJDragHydroSSA(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*return if floating (gradient is 0)*/ @@ -2316,6 +2333,5 @@ if(domaintype!=Domain2DhorizontalEnum){basalelement->DeleteMaterials(); delete basalelement;}; }/*}}}*/ - -void AdjointHorizAnalysis::GradientJDSSA(Element* element,Vector* gradient,int control_index){/*{{{*/ +void AdjointHorizAnalysis::GradientJDSSA(Element* element,Vector* gradient,int control_interp,int control_index){/*{{{*/ /*Intermediaries*/ Index: /issm/trunk-jpl/src/c/analyses/AdjointHorizAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/AdjointHorizAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/AdjointHorizAnalysis.h (revision 25317) @@ -35,25 +35,25 @@ ElementVector* CreatePVectorSSA(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); - void GradientJBbarFS(Element* element,Vector* gradient,int control_index); - void GradientJBbarGradient(Element* element,Vector* gradient,int control_index); - void GradientJBinitial(Element* element,Vector* gradient,int control_index); - void GradientJBbarL1L2(Element* element,Vector* gradient,int control_index); - void GradientJBbarHO(Element* element,Vector* gradient,int control_index); - void GradientJBbarSSA(Element* element,Vector* gradient,int control_index); - void GradientJBFS(Element* element,Vector* gradient,int control_index); - void GradientJBGradient(Element* element,Vector* gradient,int control_index); - void GradientJBHO(Element* element,Vector* gradient,int control_index); - void GradientJBSSA(Element* element,Vector* gradient,int control_index); - void GradientJDragFS(Element* element,Vector* gradient,int control_index); - void GradientJDragGradient(Element* element,Vector* gradient,int control_index); - void GradientJDragL1L2(Element* element,Vector* gradient,int control_index); - void GradientJDragHO(Element* element,Vector* gradient,int control_index); - void GradientJDragSSA(Element* element,Vector* gradient,int control_index); - void GradientJDragHydroFS(Element* element,Vector* gradient,int control_index); - void GradientJDragHydroL1L2(Element* element,Vector* gradient,int control_index); - void GradientJDragHydroHO(Element* element,Vector* gradient,int control_index); - void GradientJDragHydroSSA(Element* element,Vector* gradient,int control_index); - void GradientJDSSA(Element* element,Vector* gradient,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); + void GradientJBbarFS(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJBbarGradient(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJBinitial(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJBbarL1L2(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJBbarHO(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJBbarSSA(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJBFS(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJBGradient(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJBHO(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJBSSA(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDragFS(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDragGradient(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDragL1L2(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDragHO(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDragSSA(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDragHydroFS(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDragHydroL1L2(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDragHydroHO(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDragHydroSSA(Element* element,Vector* gradient,int control_interp,int control_index); + void GradientJDSSA(Element* element,Vector* gradient,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void InputUpdateFromSolutionFS(IssmDouble* solution,Element* element); Index: /issm/trunk-jpl/src/c/analyses/Analysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/Analysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/Analysis.h (revision 25317) @@ -49,5 +49,5 @@ virtual ElementVector* CreatePVector(Element* element)=0; virtual void GetSolutionFromInputs(Vector* solution,Element* element)=0; - virtual void GradientJ(Vector* gradient,Element* element,int control_type,int control_index)=0; + virtual void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index)=0; virtual void InputUpdateFromSolution(IssmDouble* solution,Element* element)=0; virtual void UpdateConstraints(FemModel* femmodel)=0; Index: /issm/trunk-jpl/src/c/analyses/Balancethickness2Analysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/Balancethickness2Analysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/Balancethickness2Analysis.cpp (revision 25317) @@ -174,5 +174,5 @@ element->GetSolutionFromInputsOneDof(solution,ThicknessEnum); }/*}}}*/ -void Balancethickness2Analysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void Balancethickness2Analysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/Balancethickness2Analysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/Balancethickness2Analysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/Balancethickness2Analysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/BalancethicknessAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/BalancethicknessAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/BalancethicknessAnalysis.cpp (revision 25317) @@ -424,5 +424,5 @@ _error_("not implemented yet"); }/*}}}*/ -void BalancethicknessAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void BalancethicknessAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ /* WARNING: this gradient is valid for Soft balance thickness only */ Index: /issm/trunk-jpl/src/c/analyses/BalancethicknessAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/BalancethicknessAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/BalancethicknessAnalysis.h (revision 25317) @@ -31,5 +31,5 @@ ElementVector* CreatePVectorDG(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/BalancethicknessSoftAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/BalancethicknessSoftAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/BalancethicknessSoftAnalysis.cpp (revision 25317) @@ -45,5 +45,5 @@ _error_("not implemented yet"); }/*}}}*/ -void BalancethicknessSoftAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void BalancethicknessSoftAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/BalancethicknessSoftAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/BalancethicknessSoftAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/BalancethicknessSoftAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/BalancevelocityAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/BalancevelocityAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/BalancevelocityAnalysis.cpp (revision 25317) @@ -239,5 +239,5 @@ _error_("not implemented yet"); }/*}}}*/ -void BalancevelocityAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void BalancevelocityAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/BalancevelocityAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/BalancevelocityAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/BalancevelocityAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.cpp (revision 25317) @@ -677,5 +677,5 @@ element->GetSolutionFromInputsOneDof(solution,DamageDbarEnum); }/*}}}*/ -void DamageEvolutionAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void DamageEvolutionAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.h (revision 25317) @@ -31,5 +31,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/DepthAverageAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/DepthAverageAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/DepthAverageAnalysis.cpp (revision 25317) @@ -140,5 +140,5 @@ _error_("not implemented yet"); }/*}}}*/ -void DepthAverageAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void DepthAverageAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/DepthAverageAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/DepthAverageAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/DepthAverageAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp (revision 25317) @@ -1521,5 +1521,5 @@ return (1.-waterfraction)*kappa_i + waterfraction*kappa_w; }/*}}}*/ -void EnthalpyAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void EnthalpyAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.h (revision 25317) @@ -48,5 +48,5 @@ static int GetThermalBasalCondition(Element* element, IssmDouble enthalpy, IssmDouble enthalpy_up, IssmDouble pressure, IssmDouble pressure_up, IssmDouble watercolumn, IssmDouble meltingrate); static IssmDouble GetWetIceConductivity(Element* element, IssmDouble enthalpy, IssmDouble pressure); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); static void PostProcessing(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/EnumToAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/EnumToAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/EnumToAnalysis.h (revision 25317) @@ -22,5 +22,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/EsaAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/EsaAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/EsaAnalysis.cpp (revision 25317) @@ -205,5 +205,5 @@ _error_("not implemented yet"); }/*}}}*/ -void EsaAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void EsaAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/EsaAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/EsaAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/EsaAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/ExtrapolationAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/ExtrapolationAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/ExtrapolationAnalysis.cpp (revision 25317) @@ -211,5 +211,5 @@ _error_("not implemented yet"); }/*}}}*/ -void ExtrapolationAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void ExtrapolationAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/ExtrapolationAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/ExtrapolationAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/ExtrapolationAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); int GetExtrapolationCase(Element* element); Index: /issm/trunk-jpl/src/c/analyses/ExtrudeFromBaseAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/ExtrudeFromBaseAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/ExtrudeFromBaseAnalysis.cpp (revision 25317) @@ -95,5 +95,5 @@ _error_("not implemented yet"); }/*}}}*/ -void ExtrudeFromBaseAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void ExtrudeFromBaseAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/ExtrudeFromBaseAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/ExtrudeFromBaseAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/ExtrudeFromBaseAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/ExtrudeFromTopAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/ExtrudeFromTopAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/ExtrudeFromTopAnalysis.cpp (revision 25317) @@ -95,5 +95,5 @@ _error_("not implemented yet"); }/*}}}*/ -void ExtrudeFromTopAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void ExtrudeFromTopAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/ExtrudeFromTopAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/ExtrudeFromTopAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/ExtrudeFromTopAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/FreeSurfaceBaseAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/FreeSurfaceBaseAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/FreeSurfaceBaseAnalysis.cpp (revision 25317) @@ -348,5 +348,5 @@ _error_("not implemented yet"); }/*}}}*/ -void FreeSurfaceBaseAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void FreeSurfaceBaseAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/FreeSurfaceBaseAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/FreeSurfaceBaseAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/FreeSurfaceBaseAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/FreeSurfaceTopAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/FreeSurfaceTopAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/FreeSurfaceTopAnalysis.cpp (revision 25317) @@ -322,5 +322,5 @@ _error_("not implemented yet"); }/*}}}*/ -void FreeSurfaceTopAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void FreeSurfaceTopAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/FreeSurfaceTopAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/FreeSurfaceTopAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/FreeSurfaceTopAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/GLheightadvectionAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/GLheightadvectionAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/GLheightadvectionAnalysis.cpp (revision 25317) @@ -186,5 +186,5 @@ _error_("not implemented yet"); }/*}}}*/ -void GLheightadvectionAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void GLheightadvectionAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/GLheightadvectionAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/GLheightadvectionAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/GLheightadvectionAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/GiaAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/GiaAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/GiaAnalysis.cpp (revision 25317) @@ -106,5 +106,5 @@ _error_("not implemented yet"); }/*}}}*/ -void GiaAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void GiaAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/GiaAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/GiaAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/GiaAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.cpp (revision 25317) @@ -384,5 +384,5 @@ element->GetSolutionFromInputsOneDof(solution,EplHeadSubstepEnum); }/*}}}*/ -void HydrologyDCEfficientAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void HydrologyDCEfficientAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyDCEfficientAnalysis.h (revision 25317) @@ -30,5 +30,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.cpp (revision 25317) @@ -450,5 +450,5 @@ element->GetSolutionFromInputsOneDof(solution,SedimentHeadSubstepEnum); }/*}}}*/ -void HydrologyDCInefficientAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void HydrologyDCInefficientAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyDCInefficientAnalysis.h (revision 25317) @@ -28,5 +28,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/HydrologyGlaDSAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyGlaDSAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyGlaDSAnalysis.cpp (revision 25317) @@ -430,5 +430,5 @@ }/*}}}*/ -void HydrologyGlaDSAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void HydrologyGlaDSAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/HydrologyGlaDSAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyGlaDSAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyGlaDSAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/HydrologyPismAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyPismAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyPismAnalysis.cpp (revision 25317) @@ -79,5 +79,5 @@ _error_("not implemented"); }/*}}}*/ -void HydrologyPismAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void HydrologyPismAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/HydrologyPismAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyPismAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyPismAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.cpp (revision 25317) @@ -343,5 +343,5 @@ element->GetSolutionFromInputsOneDof(solution,HydrologyHeadEnum); }/*}}}*/ -void HydrologyShaktiAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void HydrologyShaktiAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyShaktiAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/HydrologyShreveAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyShreveAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyShreveAnalysis.cpp (revision 25317) @@ -281,5 +281,5 @@ element->GetSolutionFromInputsOneDof(solution,WatercolumnEnum); }/*}}}*/ -void HydrologyShreveAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void HydrologyShreveAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/HydrologyShreveAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/HydrologyShreveAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/HydrologyShreveAnalysis.h (revision 25317) @@ -28,5 +28,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/L2ProjectionBaseAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/L2ProjectionBaseAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/L2ProjectionBaseAnalysis.cpp (revision 25317) @@ -198,5 +198,5 @@ _error_("not implemented yet"); }/*}}}*/ -void L2ProjectionBaseAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void L2ProjectionBaseAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/L2ProjectionBaseAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/L2ProjectionBaseAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/L2ProjectionBaseAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/L2ProjectionEPLAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/L2ProjectionEPLAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/L2ProjectionEPLAnalysis.cpp (revision 25317) @@ -223,5 +223,5 @@ _error_("not implemented yet"); }/*}}}*/ -void L2ProjectionEPLAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void L2ProjectionEPLAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/L2ProjectionEPLAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/L2ProjectionEPLAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/L2ProjectionEPLAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp (revision 25317) @@ -621,5 +621,5 @@ _error_("not implemented yet"); }/*}}}*/ -void LevelsetAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void LevelsetAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.h (revision 25317) @@ -28,5 +28,5 @@ IssmDouble GetDistanceToStraight(IssmDouble* q, IssmDouble* s0, IssmDouble* s1); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/LoveAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/LoveAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/LoveAnalysis.cpp (revision 25317) @@ -67,5 +67,5 @@ _error_("not implemented yet"); }/*}}}*/ -void LoveAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void LoveAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/LoveAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/LoveAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/LoveAnalysis.h (revision 25317) @@ -29,5 +29,5 @@ void GetBprime(IssmDouble* B,Element* element,int dim,IssmDouble* xyz_list,Gauss* gauss); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp (revision 25317) @@ -790,5 +790,5 @@ element->GetSolutionFromInputsOneDof(solution,ThicknessEnum); }/*}}}*/ -void MasstransportAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void MasstransportAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.h (revision 25317) @@ -31,5 +31,5 @@ ElementVector* CreatePVectorDG(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/MeltingAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/MeltingAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/MeltingAnalysis.cpp (revision 25317) @@ -125,5 +125,5 @@ _error_("not implemented yet"); }/*}}}*/ -void MeltingAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void MeltingAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/MeltingAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/MeltingAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/MeltingAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/SealevelriseAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/SealevelriseAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/SealevelriseAnalysis.cpp (revision 25317) @@ -413,5 +413,5 @@ _error_("not implemented yet"); }/*}}}*/ -void SealevelriseAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void SealevelriseAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/SealevelriseAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/SealevelriseAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/SealevelriseAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/SmbAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/SmbAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/SmbAnalysis.cpp (revision 25317) @@ -507,5 +507,5 @@ _error_("not implemented yet"); }/*}}}*/ -void SmbAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void SmbAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/SmbAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/SmbAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/SmbAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/SmoothAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/SmoothAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/SmoothAnalysis.cpp (revision 25317) @@ -206,5 +206,5 @@ _error_("not implemented yet"); }/*}}}*/ -void SmoothAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void SmoothAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/SmoothAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/SmoothAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/SmoothAnalysis.h (revision 25317) @@ -27,5 +27,5 @@ ElementVector* CreatePVector(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.cpp (revision 25317) @@ -1237,5 +1237,5 @@ xDelete(doflist); }/*}}}*/ -void StressbalanceAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void StressbalanceAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/StressbalanceAnalysis.h (revision 25317) @@ -28,5 +28,5 @@ void GetSolutionFromInputs(Vector* solution,Element* element); void GetSolutionFromInputsHoriz(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/StressbalanceSIAAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/StressbalanceSIAAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/StressbalanceSIAAnalysis.cpp (revision 25317) @@ -556,5 +556,5 @@ xDelete(values); }/*}}}*/ -void StressbalanceSIAAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void StressbalanceSIAAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/StressbalanceSIAAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/StressbalanceSIAAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/StressbalanceSIAAnalysis.h (revision 25317) @@ -31,5 +31,5 @@ ElementVector* CreatePVector3D(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/StressbalanceVerticalAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/StressbalanceVerticalAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/StressbalanceVerticalAnalysis.cpp (revision 25317) @@ -505,5 +505,5 @@ element->GetSolutionFromInputsOneDof(solution,VzEnum); }/*}}}*/ -void StressbalanceVerticalAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void StressbalanceVerticalAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/StressbalanceVerticalAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/StressbalanceVerticalAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/StressbalanceVerticalAnalysis.h (revision 25317) @@ -33,5 +33,5 @@ ElementVector* CreatePVectorVolume(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/ThermalAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/ThermalAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/ThermalAnalysis.cpp (revision 25317) @@ -742,5 +742,5 @@ element->GetSolutionFromInputsOneDof(solution,TemperatureEnum); }/*}}}*/ -void ThermalAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void ThermalAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/ThermalAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/ThermalAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/ThermalAnalysis.h (revision 25317) @@ -32,5 +32,5 @@ ElementVector* CreatePVectorVolume(Element* element); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/analyses/UzawaPressureAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/UzawaPressureAnalysis.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/UzawaPressureAnalysis.cpp (revision 25317) @@ -182,5 +182,5 @@ _error_("not implemented yet"); }/*}}}*/ -void UzawaPressureAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_index){/*{{{*/ +void UzawaPressureAnalysis::GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index){/*{{{*/ _error_("Not implemented yet"); }/*}}}*/ Index: /issm/trunk-jpl/src/c/analyses/UzawaPressureAnalysis.h =================================================================== --- /issm/trunk-jpl/src/c/analyses/UzawaPressureAnalysis.h (revision 25316) +++ /issm/trunk-jpl/src/c/analyses/UzawaPressureAnalysis.h (revision 25317) @@ -28,5 +28,5 @@ void GetM(IssmDouble* M,Element* element,Gauss* gauss); void GetSolutionFromInputs(Vector* solution,Element* element); - void GradientJ(Vector* gradient,Element* element,int control_type,int control_index); + void GradientJ(Vector* gradient,Element* element,int control_type,int control_interp,int control_index); void InputUpdateFromSolution(IssmDouble* solution,Element* element); void UpdateConstraints(FemModel* femmodel); Index: /issm/trunk-jpl/src/c/classes/Elements/Element.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Element.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Elements/Element.cpp (revision 25317) @@ -1509,35 +1509,39 @@ int num_controls; bool isautodiff; + int* N=NULL; + int* M=NULL; + int* interp=NULL; parameters->FindParam(&num_controls,InversionNumControlParametersEnum); - parameters->FindParam(&isautodiff,AutodiffIsautodiffEnum); + parameters->FindParam(&N,NULL,ControlInputSizeNEnum); + parameters->FindParam(&M,NULL,ControlInputSizeMEnum); + parameters->FindParam(&interp,NULL,ControlInputInterpolationEnum); /*Get number of vertices*/ const int NUM_VERTICES = this->GetNumberOfVertices(); - if(isautodiff){ - int* N=NULL; - int* M=NULL; - int start = 0; - parameters->FindParam(&N,NULL,ControlInputSizeNEnum); - parameters->FindParam(&M,NULL,ControlInputSizeMEnum); - if(control_index>0) { - for(int n=0;nvertices[i]->Sid() + start + n*M[control_index]; + indexing[i+n*NUM_VERTICES]= start + n*M[control_index] + this->vertices[i]->Sid(); } } } + else if(interp[control_index]==P0Enum){ + for(int n=0;nSid(); + } + } else{ - int M; - parameters->FindParam(&M,ControlInputSizeMEnum); - /*get gradient indices*/ - for(int i=0;ivertices[i]->Sid() + (control_index)*M; - } - } + _error_("interpolation not supported"); + } + + xDelete(M); + xDelete(N); + xDelete(interp); } /*}}}*/ @@ -1776,4 +1780,10 @@ } this->AddControlInput(input_enum,inputs2,iomodel,values,values_min,values_max,P1Enum,id); + } + else if(M==iomodel->numberofelements && N==1){ + values[0] = vector[this->Sid()]; + values_min[0] = scale*min_vector[this->Sid()]; + values_max[0] = scale*max_vector[this->Sid()]; + this->AddControlInput(input_enum,inputs2,iomodel,values,values_min,values_max,P0Enum,id); } Index: /issm/trunk-jpl/src/c/classes/Elements/Element.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Element.h (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Elements/Element.h (revision 25317) @@ -231,7 +231,7 @@ virtual void ComputeEsaStrainAndVorticity(void)=0; virtual void Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters,Inputs2* inputs2in)=0; - virtual void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int N,int M)=0; + virtual void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int M,int N,int interp)=0; virtual void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index)=0; - virtual void ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum)=0; + virtual void ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum,int control_interp)=0; virtual void CreateDistanceInputFromSegmentlist(IssmDouble* distances,int distanceenum){_error_("not implemented yet");}; virtual void CreateInputTimeAverage(int transientinput_enum,int averagedinput_enum,IssmDouble init_time,IssmDouble end_time,int averaging_method){_error_("not implemented yet "<ObjectEnum());}; @@ -335,5 +335,5 @@ virtual void SetElementInput(Inputs2* inputs2,int enum_in,IssmDouble values){_error_("not implemented yet");}; virtual void SetElementInput(Inputs2* inputs2,int numindices,int* indices,IssmDouble* values,int enum_in){_error_("not implemented yet");}; - virtual void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset,int N,int M)=0; + virtual void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset,int M,int N)=0; virtual void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index)=0; virtual void SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Materials* materials,Parameters* parameters)=0; Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 25317) @@ -899,33 +899,38 @@ } /*}}}*/ -void Penta::ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int N,int M){/*{{{*/ - - if(enum_type==MaterialsRheologyBbarEnum) enum_type = MaterialsRheologyBEnum; - if(enum_type==DamageDbarEnum) enum_type = DamageDEnum; - - _error_("not implemented"); - int vertexpidlist[NUMVERTICES]; - IssmDouble grad_list[NUMVERTICES]; - Input2* grad_input=NULL; - Input2* input=NULL; - - if(enum_type==MaterialsRheologyBbarEnum){ - input=this->GetInput2(MaterialsRheologyBEnum); - } - else if(enum_type==DamageDbarEnum){ - input=this->GetInput2(DamageDEnum); - } - else{ - input=this->GetInput2(enum_type); - } - if (!input) _error_("Input " << EnumToStringx(enum_type) << " not found"); - if(input->ObjectEnum()!=ControlInput2Enum) _error_("Input " << EnumToStringx(enum_type) << " is not a ControlInput"); - - GradientIndexing(&vertexpidlist[0],control_index); - - //for(int i=0;iSetGradient(grad_input); - _error_("not implemented"); +void Penta::ControlInputSetGradient(IssmDouble* gradient,int control_enum,int control_index,int offset,int M,int N,int interp){/*{{{*/ + + IssmDouble values[NUMVERTICES]; + int lidlist[NUMVERTICES]; + int idlist[NUMVERTICES]; + + if(control_enum==MaterialsRheologyBbarEnum) control_enum = MaterialsRheologyBEnum; + if(control_enum==DamageDbarEnum) control_enum = DamageDEnum; + + ElementInput2* input=this->inputs2->GetControlInput2Data(control_enum,"gradient"); _assert_(input); + this->GetVerticesLidList(&lidlist[0]); + GradientIndexing(&idlist[0],control_index); + + /*Get values on vertices*/ + if(input->ObjectEnum()==PentaInput2Enum && input->GetInputInterpolationType()==P1Enum){ + _assert_(N==1); + for(int i=0;iSetInput(P1Enum,NUMVERTICES,&lidlist[0],&values[0]); + } + else if(input->ObjectEnum()==PentaInput2Enum && input->GetInputInterpolationType()==P0Enum){ + _assert_(N==1); + input->SetInput(P0Enum,this->lid,gradient[idlist[0]]); + } + else if(input->ObjectEnum()==TransientInput2Enum){ + for(int n=0;nSetInput(new_input,n); + //controlinput->Configure(parameters); + } + } + else _error_("Type not supported"); } @@ -946,5 +951,5 @@ this->inputs2->SetTriaControlInputGradient(enum_type,P1Enum,NUMVERTICES,&vertexlids[0],&grad_list[0]); }/*}}}*/ -void Penta::ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum){/*{{{*/ +void Penta::ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum,int control_interp){/*{{{*/ int sidlist[NUMVERTICES]; @@ -955,30 +960,50 @@ IssmDouble value,gradient; + /*Get relevant inputs*/ if(control_enum==MaterialsRheologyBbarEnum) control_enum = MaterialsRheologyBEnum; if(control_enum==DamageDbarEnum) control_enum = DamageDEnum; - - this->GetVerticesConnectivityList(&connectivity[0]); - this->GetVerticesSidList(&sidlist[0]); - this->GetVerticesLidList(&lidlist[0]); - ElementInput2* control_value = this->inputs2->GetControlInput2Data(control_enum,"value"); _assert_(control_value); ElementInput2* control_gradient = this->inputs2->GetControlInput2Data(control_enum,"gradient"); _assert_(control_gradient); - control_value->Serve(NUMVERTICES,&lidlist[0]); - control_gradient->Serve(NUMVERTICES,&lidlist[0]); - - GaussPenta* gauss=new GaussPenta(); - for (int iv=0;ivGaussVertex(iv); - - control_value->GetInputValue(&value,gauss); - control_gradient->GetInputValue(&gradient,gauss); - - values[iv] = reCast(value)/reCast(connectivity[iv]); - gradients[iv] = reCast(gradient)/reCast(connectivity[iv]); - } - delete gauss; - - vector_control->SetValues(NUMVERTICES,&sidlist[0],&values[0],ADD_VAL); - vector_gradient->SetValues(NUMVERTICES,&sidlist[0],&gradients[0],ADD_VAL); + + if(control_interp==P1Enum){ + _assert_(control_value->GetInputInterpolationType()==P1Enum); + _assert_(control_gradient->GetInputInterpolationType()==P1Enum); + + this->GetVerticesConnectivityList(&connectivity[0]); + this->GetVerticesSidList(&sidlist[0]); + this->GetVerticesLidList(&lidlist[0]); + + control_value->Serve(NUMVERTICES,&lidlist[0]); + control_gradient->Serve(NUMVERTICES,&lidlist[0]); + + GaussPenta* gauss=new GaussPenta(); + for (int iv=0;ivGaussVertex(iv); + + control_value->GetInputValue(&value,gauss); + control_gradient->GetInputValue(&gradient,gauss); + + values[iv] = reCast(value)/reCast(connectivity[iv]); + gradients[iv] = reCast(gradient)/reCast(connectivity[iv]); + } + delete gauss; + + vector_control->SetValues(NUMVERTICES,&sidlist[0],&values[0],ADD_VAL); + vector_gradient->SetValues(NUMVERTICES,&sidlist[0],&gradients[0],ADD_VAL); + } + else if(control_interp==P0Enum){ + _assert_(control_value->GetInputInterpolationType()==P0Enum); + _assert_(control_gradient->GetInputInterpolationType()==P0Enum); + + control_value->Serve(1,&this->lid); + control_gradient->Serve(1,&this->lid); + + vector_control->SetValue(this->sid,control_value->element_values[0],ADD_VAL); + vector_gradient->SetValue(this->sid,control_gradient->element_values[0],ADD_VAL); + } + else{ + _error_("not supported"); + } + }/*}}}*/ void Penta::CreateDistanceInputFromSegmentlist(IssmDouble* distances,int distanceenum){/*{{{*/ @@ -1822,80 +1847,136 @@ void Penta::GetVectorFromControlInputs(Vector* vector,int control_enum,int control_index,const char* data){/*{{{*/ - /*Get out if this is not an element input*/ - if(!IsInputEnum(control_enum)) _error_("Enum "<inputs2->GetControlInput2Data(control_enum,data); _assert_(input); - /*Intermediaries*/ - int numindices; - int indices[NUMVERTICES]; - - /*Check interpolation*/ - int interpolation = input->GetInterpolation(); - switch(interpolation){ - case P1Enum: - numindices = NUMVERTICES; - for(int i=0;ilid; - input->Serve(numindices,&indices[0]); - break; - default: _error_("interpolation "<NewGauss(); - for(int iv=0;ivGaussVertex(iv); - input->GetInputValue(&values[iv],gauss); - } - delete gauss; - - vector->SetValues(NUMVERTICES,idlist,&values[0],INS_VAL); -} -/*}}}*/ -void Penta::GetVectorFromControlInputs(Vector* vector,int control_enum,int control_index,const char* data,int offset){/*{{{*/ - - int* idlist = NULL; - IssmDouble* values = NULL; - - /*Get input*/ - ElementInput2* input=this->inputs2->GetControlInput2Data(control_enum,data); _assert_(input); + /*Lid list once for all*/ + int lidlist[NUMVERTICES]; + for(int i=0;ilid; /*Check what input we are dealing with*/ switch(input->ObjectEnum()){ case PentaInput2Enum: - { - PentaInput2* pentainput = xDynamicCast(input); - if(pentainput->GetInputInterpolationType()!=P1Enum) _error_("not supported yet"); - - /*Create list of indices and values for global vector*/ - idlist = xNew(NUMVERTICES); - values = xNew(NUMVERTICES); - GradientIndexing(&idlist[0],control_index); - for(int i=0;ielement_values[i]; + { + IssmDouble values[NUMVERTICES]; + int idlist[NUMVERTICES]; + + PentaInput2* triainput = xDynamicCast(input); + + /*Create list of indices and values for global vector*/ + GradientIndexing(&idlist[0],control_index); + + if(triainput->GetInputInterpolationType()==P1Enum){ + input->Serve(NUMVERTICES,&lidlist[0]); + for(int i=0;ielement_values[i]; + vector->SetValues(NUMVERTICES,idlist,values,INS_VAL); + } + else if(triainput->GetInputInterpolationType()==P0Enum){ + input->Serve(1,&this->lid); + vector->SetValue(idlist[0],triainput->element_values[0],INS_VAL); + } + else{ + _error_("not supported yet"); + } + break; + } + + case TransientInputEnum: + { + TransientInput2* transientinput = xDynamicCast(input); + int N = transientinput->numtimesteps; + int* M = NULL; + parameters->FindParam(&M,NULL,ControlInputSizeMEnum); + int* idlist = xNew(NUMVERTICES*N); + IssmDouble* values = xNew(NUMVERTICES*N); + for(int t=0;tnumtimesteps;t++) { + IssmDouble time = transientinput->GetTimeByOffset(t); + _error_("not implemented"); + //PentaInput* timeinput = xDynamicCast(transientinput->GetTimeInput(time)); + //if(timeinput->interpolation_type!=P1Enum) _error_("not supported yet"); + //input->Serve(NUMVERTICES,&lidlist[0]); + ///*Create list of indices and values for global vector*/ + //for(int i=0;ivertices[i]->Sid()+t*M[control_index]; + // values[N*i+t] = timeinput->values[i]; + //} } + vector->SetValues(NUMVERTICES*transientinput->numtimesteps,idlist,values,INS_VAL); + xDelete(M); + xDelete(idlist); + xDelete(values); + break; + } + default: _error_("input "<ObjectEnum())<<" not supported yet"); + } +} +/*}}}*/ +void Penta::GetVectorFromControlInputs(Vector* vector,int control_enum,int control_index,const char* data,int offset){/*{{{*/ + + /*Get input*/ + if(control_enum==MaterialsRheologyBbarEnum) control_enum=MaterialsRheologyBEnum; + ElementInput2* input=this->inputs2->GetControlInput2Data(control_enum,data); _assert_(input); + + /*Lid list once for all*/ + int lidlist[NUMVERTICES]; + for(int i=0;ilid; + + /*Check what input we are dealing with*/ + switch(input->ObjectEnum()){ + case PentaInput2Enum: + { + IssmDouble values[NUMVERTICES]; + int idlist[NUMVERTICES]; + + PentaInput2* triainput = xDynamicCast(input); + + /*Create list of indices and values for global vector*/ + GradientIndexing(&idlist[0],control_index); + + if(triainput->GetInputInterpolationType()==P1Enum){ + input->Serve(NUMVERTICES,&lidlist[0]); + for(int i=0;ielement_values[i]; vector->SetValues(NUMVERTICES,idlist,values,INS_VAL); + } + else if(triainput->GetInputInterpolationType()==P0Enum){ + input->Serve(1,&this->lid); + vector->SetValue(idlist[0],triainput->element_values[0],INS_VAL); + } + else{ + _error_("not supported yet"); + } + break; + } + + case TransientInputEnum: + { + TransientInput2* transientinput = xDynamicCast(input); + int N = transientinput->numtimesteps; + int* M = NULL; + parameters->FindParam(&M,NULL,ControlInputSizeMEnum); + int* idlist = xNew(NUMVERTICES*N); + IssmDouble* values = xNew(NUMVERTICES*N); + for(int t=0;tnumtimesteps;t++) { + IssmDouble time = transientinput->GetTimeByOffset(t); + _error_("not implemented"); + //PentaInput* timeinput = xDynamicCast(transientinput->GetTimeInput(time)); + //if(timeinput->interpolation_type!=P1Enum) _error_("not supported yet"); + //input->Serve(NUMVERTICES,&lidlist[0]); + ///*Create list of indices and values for global vector*/ + //for(int i=0;ivertices[i]->Sid()+t*M[control_index]; + // values[N*i+t] = timeinput->values[i]; + //} + } + vector->SetValues(NUMVERTICES*transientinput->numtimesteps,idlist,values,INS_VAL); + xDelete(M); + xDelete(idlist); + xDelete(values); break; - } - - case TransientInputEnum: - { - _error_("not implemented (see Tria)"); - break; - } - default: _error_("input "<ObjectEnum()<<" not supported yet"); - - } - xDelete(idlist); - xDelete(values); - } -/*}}}*/ + } + default: _error_("input "<ObjectEnum())<<" not supported yet"); + } +}/*}}}*/ void Penta::GetVerticesCoordinatesBase(IssmDouble** pxyz_list){/*{{{*/ @@ -3303,8 +3384,9 @@ } /*}}}*/ -void Penta::SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset, int N, int M){/*{{{*/ +void Penta::SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset, int M, int N){/*{{{*/ IssmDouble values[NUMVERTICES]; - int vertexpidlist[NUMVERTICES],control_init; + int lidlist[NUMVERTICES]; + int idlist[NUMVERTICES],control_init; /*Specific case for depth averaged quantities*/ @@ -3319,32 +3401,62 @@ } + /*Get Domain type*/ + int domaintype; + parameters->FindParam(&domaintype,DomainTypeEnum); + + /*Specific case for depth averaged quantities*/ + if(domaintype==Domain2DverticalEnum){ + if(control_enum==MaterialsRheologyBbarEnum){ + control_enum=MaterialsRheologyBEnum; + if(!IsOnBase()) return; + } + if(control_enum==DamageDbarEnum){ + control_enum=DamageDEnum; + if(!IsOnBase()) return; + } + } + /*Get out if this is not an element input*/ if(!IsInputEnum(control_enum)) return; /*Prepare index list*/ - GradientIndexing(&vertexpidlist[0],control_index); + ElementInput2* input=this->inputs2->GetControlInput2Data(control_enum,"value"); _assert_(input); + this->GetVerticesLidList(&lidlist[0]); + GradientIndexing(&idlist[0],control_index); /*Get values on vertices*/ - for(int i=0;iGetInput2(control_enum); _assert_(input); - if(input->ObjectEnum()!=ControlInputEnum){ - _error_("input " << EnumToStringx(control_enum) << " is not a ControlInput"); - } - - //((ControlInput*)input)->SetInput(new_input); - + if(input->ObjectEnum()==PentaInput2Enum && input->GetInputInterpolationType()==P1Enum){ + _assert_(N==1); + for(int i=0;iSetInput(P1Enum,NUMVERTICES,&lidlist[0],&values[0]); + } + else if(input->ObjectEnum()==PentaInput2Enum && input->GetInputInterpolationType()==P0Enum){ + _assert_(N==1); + input->SetInput(P0Enum,this->lid,vector[idlist[0]]); + } + else if(input->ObjectEnum()==TransientInput2Enum){ + for(int n=0;nSetInput(new_input,n); + //controlinput->Configure(parameters); + } + } + else _error_("Type not supported"); + + /*Extrude depending on the control*/ if(control_init==MaterialsRheologyBbarEnum){ - this->InputExtrude(control_enum,-1); + this->ControlInputExtrude(control_enum,-1); } if(control_init==DamageDbarEnum){ - this->InputExtrude(control_enum,-1); + this->ControlInputExtrude(control_enum,-1); } } /*}}}*/ void Penta::SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index){/*{{{*/ + + _error_("not needed anymore?"); IssmDouble values[NUMVERTICES]; Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Penta.h (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Elements/Penta.h (revision 25317) @@ -63,7 +63,7 @@ void ComputeStressTensor(); void Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters,Inputs2* inputs2in); - void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int N,int M); + void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int M,int N,int interp); void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index); - void ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum); + void ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum,int control_interp); void CreateDistanceInputFromSegmentlist(IssmDouble* distances,int distanceenum); ElementMatrix* CreateBasalMassMatrix(void); @@ -167,5 +167,5 @@ void SetElementInput(Inputs2* inputs2,int enum_in,IssmDouble values); void SetElementInput(Inputs2* inputs2,int numindices,int* indices,IssmDouble* values,int enum_in); - void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset, int N,int M); + void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset, int M,int N); void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index); void SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Materials* materials,Parameters* parameters); Index: /issm/trunk-jpl/src/c/classes/Elements/Seg.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Seg.h (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Elements/Seg.h (revision 25317) @@ -50,7 +50,7 @@ void ComputeStressTensor(){_error_("not implemented yet");}; void Configure(Elements* elements,Loads* loads,Nodes* nodesin,Vertices* verticesin,Materials* materials,Parameters* parameters,Inputs2* inputs2in){_error_("not implemented yet");}; - void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int N,int M){_error_("not implemented yet");}; + void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int M,int N,int interp){_error_("not implemented yet");}; void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index){_error_("not implemented yet");}; - void ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum){_error_("not implemented yet");}; + void ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum,int control_interp){_error_("not implemented yet");}; void ElementResponse(IssmDouble* presponse,int response_enum){_error_("not implemented yet");}; void ElementSizes(IssmDouble* hx,IssmDouble* hy,IssmDouble* hz){_error_("not implemented yet");}; @@ -134,5 +134,5 @@ void ResetFSBasalBoundaryCondition(void){_error_("not implemented yet");}; void ResetHooks(){_error_("not implemented yet");}; - void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset,int N,int M){_error_("not implemented yet");}; + void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset,int M,int N){_error_("not implemented yet");}; void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index){_error_("not implemented yet");}; void SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Materials* materials,Parameters* parameters){_error_("not implemented yet");}; Index: /issm/trunk-jpl/src/c/classes/Elements/Tetra.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tetra.h (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Elements/Tetra.h (revision 25317) @@ -48,7 +48,7 @@ void ComputeEsaStrainAndVorticity(){_error_("not implemented yet!");}; void Configure(Elements* elements,Loads* loads,Nodes* nodesin,Vertices* verticesin,Materials* materials,Parameters* parameters,Inputs2* inputs2in); - void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int N,int M){_error_("not implemented yet");}; + void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int M,int N,int interp){_error_("not implemented yet");}; void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index){_error_("not implemented yet");}; - void ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum){_error_("not implemented yet");}; + void ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum,int control_interp){_error_("not implemented yet");}; IssmDouble DragCoefficientAbsGradient(void){_error_("not implemented yet");}; void ElementResponse(IssmDouble* presponse,int response_enum){_error_("not implemented yet");}; @@ -140,5 +140,5 @@ void ResetHooks(); void ReduceMatrices(ElementMatrix* Ke,ElementVector* pe); - void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset,int N, int M){_error_("not implemented yet");}; + void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset,int M, int N){_error_("not implemented yet");}; void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index){_error_("not implemented yet");}; void SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Materials* materials,Parameters* parameters); Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 25317) @@ -222,4 +222,7 @@ /*Call inputs method*/ switch(interpolation_enum){ + case P0Enum: + inputs2->SetTriaControlInput(input_enum,TriaInput2Enum,interpolation_enum,id,1,&this->lid,values,values_min,values_max); + break; case P1Enum: inputs2->SetTriaControlInput(input_enum,TriaInput2Enum,interpolation_enum,id,NUMVERTICES,vertexlids,values,values_min,values_max); @@ -1097,29 +1100,39 @@ this->inputs2=inputs2in; }/*}}}*/ -void Tria::ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int N, int M){/*{{{*/ - +void Tria::ControlInputSetGradient(IssmDouble* gradient,int control_enum,int control_index,int offset,int M,int N,int interp){/*{{{*/ + + IssmDouble values[NUMVERTICES]; + int lidlist[NUMVERTICES]; int idlist[NUMVERTICES]; - int vertexlids[NUMVERTICES]; - int gradidlist[NUMVERTICES]; - IssmDouble grad_list[NUMVERTICES]; - - GradientIndexing(&gradidlist[0],control_index); - for(int i=0;ivertices[i]->lid; - - if(N==1){ - this->inputs2->SetTriaControlInputGradient(enum_type,P1Enum,NUMVERTICES,&vertexlids[0],&grad_list[0]); - } - else{ + + ElementInput2* input=this->inputs2->GetControlInput2Data(control_enum,"gradient"); _assert_(input); + this->GetVerticesLidList(&lidlist[0]); + GradientIndexing(&idlist[0],control_index); + + /*Get values on vertices*/ + if(input->ObjectEnum()==TriaInput2Enum && input->GetInputInterpolationType()==P1Enum){ + _assert_(N==1); + for(int i=0;iSetInput(P1Enum,NUMVERTICES,&lidlist[0],&values[0]); + } + else if(input->ObjectEnum()==TriaInput2Enum && input->GetInputInterpolationType()==P0Enum){ + _assert_(N==1); + input->SetInput(P0Enum,this->lid,gradient[idlist[0]]); + } + else if(input->ObjectEnum()==TransientInput2Enum){ for(int n=0;nvertices[i]->Sid()+n*M; - grad_list[i]=gradient[idlist[i]]; - } - this->inputs2->SetTriaControlInputGradient(enum_type,P1Enum,NUMVERTICES,&vertexlids[0],&grad_list[0],n); - } - } + _error_("not implemented"); + //Input* new_input = new TriaInput(control_enum,gradient,P1Enum); + //controlinput->SetInput(new_input,n); + //controlinput->Configure(parameters); + } + } + else _error_("Type not supported"); }/*}}}*/ void Tria::ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index){/*{{{*/ + _error_("NOT NEEDED ANYMORE (remove)"); int idlist[NUMVERTICES]; @@ -1134,5 +1147,5 @@ }/*}}}*/ -void Tria::ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum){/*{{{*/ +void Tria::ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum,int control_interp){/*{{{*/ int sidlist[NUMVERTICES]; @@ -1143,27 +1156,47 @@ IssmDouble value,gradient; - this->GetVerticesConnectivityList(&connectivity[0]); - this->GetVerticesSidList(&sidlist[0]); - this->GetVerticesLidList(&lidlist[0]); - + /*Get relevant inputs*/ ElementInput2* control_value = this->inputs2->GetControlInput2Data(control_enum,"value"); _assert_(control_value); ElementInput2* control_gradient = this->inputs2->GetControlInput2Data(control_enum,"gradient"); _assert_(control_gradient); - control_value->Serve(NUMVERTICES,&lidlist[0]); - control_gradient->Serve(NUMVERTICES,&lidlist[0]); - - GaussTria* gauss=new GaussTria(); - for (int iv=0;ivGaussVertex(iv); - - control_value->GetInputValue(&value,gauss); - control_gradient->GetInputValue(&gradient,gauss); - - values[iv] = reCast(value)/reCast(connectivity[iv]); - gradients[iv] = reCast(gradient)/reCast(connectivity[iv]); - } - delete gauss; - - vector_control->SetValues(NUMVERTICES,&sidlist[0],&values[0],ADD_VAL); - vector_gradient->SetValues(NUMVERTICES,&sidlist[0],&gradients[0],ADD_VAL); + + if(control_interp==P1Enum){ + _assert_(control_value->GetInputInterpolationType()==P1Enum); + _assert_(control_gradient->GetInputInterpolationType()==P1Enum); + + this->GetVerticesConnectivityList(&connectivity[0]); + this->GetVerticesSidList(&sidlist[0]); + this->GetVerticesLidList(&lidlist[0]); + + control_value->Serve(NUMVERTICES,&lidlist[0]); + control_gradient->Serve(NUMVERTICES,&lidlist[0]); + + GaussTria* gauss=new GaussTria(); + for (int iv=0;ivGaussVertex(iv); + + control_value->GetInputValue(&value,gauss); + control_gradient->GetInputValue(&gradient,gauss); + + values[iv] = reCast(value)/reCast(connectivity[iv]); + gradients[iv] = reCast(gradient)/reCast(connectivity[iv]); + } + delete gauss; + + vector_control->SetValues(NUMVERTICES,&sidlist[0],&values[0],ADD_VAL); + vector_gradient->SetValues(NUMVERTICES,&sidlist[0],&gradients[0],ADD_VAL); + } + else if(control_interp==P0Enum){ + _assert_(control_value->GetInputInterpolationType()==P0Enum); + _assert_(control_gradient->GetInputInterpolationType()==P0Enum); + + control_value->Serve(1,&this->lid); + control_gradient->Serve(1,&this->lid); + + vector_control->SetValue(this->sid,control_value->element_values[0],ADD_VAL); + vector_gradient->SetValue(this->sid,control_gradient->element_values[0],ADD_VAL); + } + else{ + _error_("not supported"); + } }/*}}}*/ @@ -2390,4 +2423,6 @@ void Tria::GetVectorFromControlInputs(Vector* vector,int control_enum,int control_index,const char* data){/*{{{*/ + _error_("DO WE NEED THIS FUNCTION?? REALLY?? see below"); + /*Get out if this is not an element input*/ if(!IsInputEnum(control_enum)) _error_("Enum "<GetInterpolation(); switch(interpolation){ + case P0Enum: + input->Serve(1,&this->lid); + break; case P1Enum: numindices = NUMVERTICES; @@ -2449,11 +2487,20 @@ TriaInput2* triainput = xDynamicCast(input); - if(triainput->GetInputInterpolationType()!=P1Enum) _error_("not supported yet"); - input->Serve(NUMVERTICES,&lidlist[0]); /*Create list of indices and values for global vector*/ GradientIndexing(&idlist[0],control_index); - for(int i=0;ielement_values[i]; - vector->SetValues(NUMVERTICES,idlist,values,INS_VAL); + + if(triainput->GetInputInterpolationType()==P1Enum){ + input->Serve(NUMVERTICES,&lidlist[0]); + for(int i=0;ielement_values[i]; + vector->SetValues(NUMVERTICES,idlist,values,INS_VAL); + } + else if(triainput->GetInputInterpolationType()==P0Enum){ + input->Serve(1,&this->lid); + vector->SetValue(idlist[0],triainput->element_values[0],INS_VAL); + } + else{ + _error_("not supported yet"); + } break; } @@ -2462,6 +2509,6 @@ { TransientInput2* transientinput = xDynamicCast(input); - int N = transientinput->numtimesteps; - int* M=NULL; + int N = transientinput->numtimesteps; + int* M = NULL; parameters->FindParam(&M,NULL,ControlInputSizeMEnum); int* idlist = xNew(NUMVERTICES*N); @@ -3907,9 +3954,9 @@ } /*}}}*/ -void Tria::SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset,int N, int M){/*{{{*/ +void Tria::SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset,int M, int N){/*{{{*/ IssmDouble values[NUMVERTICES]; int lidlist[NUMVERTICES]; - int idlist[NUMVERTICES],control_init; + int idlist[NUMVERTICES]; /*Get Domain type*/ @@ -3918,5 +3965,5 @@ /*Specific case for depth averaged quantities*/ - control_init=control_enum; + int control_init=control_enum; if(domaintype==Domain2DverticalEnum){ if(control_enum==MaterialsRheologyBbarEnum){ @@ -3933,17 +3980,22 @@ if(!IsInputEnum(control_enum)) return; + ElementInput2* input=this->inputs2->GetControlInput2Data(control_enum,"value"); _assert_(input); this->GetVerticesLidList(&lidlist[0]); - ElementInput2* input=this->inputs2->GetControlInput2Data(control_enum,"value"); _assert_(input); + GradientIndexing(&idlist[0],control_index); /*Get values on vertices*/ - for(int n=0;nObjectEnum()==TriaInput2Enum && input->GetInputInterpolationType()==P1Enum){ + _assert_(N==1); for(int i=0;ivertices[i]->Sid()+n*M; - values[i]=vector[idlist[i]]; - } - if(input->ObjectEnum()==TriaInput2Enum){ - input->SetInput(P1Enum,NUMVERTICES,&lidlist[0],&values[0]); - } - else if(input->ObjectEnum()==TransientInput2Enum){ + values[i] = vector[idlist[i]]; + } + input->SetInput(P1Enum,NUMVERTICES,&lidlist[0],&values[0]); + } + else if(input->ObjectEnum()==TriaInput2Enum && input->GetInputInterpolationType()==P0Enum){ + _assert_(N==1); + input->SetInput(P0Enum,this->lid,vector[idlist[0]]); + } + else if(input->ObjectEnum()==TransientInput2Enum){ + for(int n=0;nConfigure(parameters); } - else _error_("Type not supported"); - } + } + else _error_("Type not supported"); } /*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tria.h (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Elements/Tria.h (revision 25317) @@ -66,7 +66,7 @@ void ComputeSurfaceNormalVelocity(); void Configure(Elements* elements,Loads* loads,Nodes* nodesin,Vertices* verticesin,Materials* materials,Parameters* parameters,Inputs2* inputs2in); - void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int N, int M); + void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index,int offset,int M,int N,int interp); void ControlInputSetGradient(IssmDouble* gradient,int enum_type,int control_index); - void ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum); + void ControlToVectors(Vector* vector_control, Vector* vector_gradient,int control_enum,int control_interp); void CreateDistanceInputFromSegmentlist(IssmDouble* distances,int distanceenum); int EdgeOnBaseIndex(); @@ -130,5 +130,5 @@ void SetElementInput(Inputs2* inputs2,int enum_in,IssmDouble values); void SetElementInput(Inputs2* inputs2,int numindices,int* indices,IssmDouble* values,int enum_in); - void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset,int N,int M); + void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index,int offset,int M,int N); void SetControlInputsFromVector(IssmDouble* vector,int control_enum,int control_index); void SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Materials* materials,Parameters* parameters); Index: /issm/trunk-jpl/src/c/classes/FemModel.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/FemModel.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/classes/FemModel.cpp (revision 25317) @@ -2174,5 +2174,8 @@ int num_controls,step; IssmDouble time; - int *control_type = NULL; + int *control_type = NULL; + int *control_interp = NULL; + int *M = NULL; + int *N = NULL; /*recover results*/ @@ -2183,4 +2186,7 @@ this->parameters->FindParam(&num_controls,InversionNumControlParametersEnum); this->parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); + this->parameters->FindParam(&M,NULL,ControlInputSizeMEnum); + this->parameters->FindParam(&N,NULL,ControlInputSizeNEnum); + this->parameters->FindParam(&control_interp,NULL,ControlInputInterpolationEnum); this->parameters->FindParam(&step,StepEnum); this->parameters->FindParam(&time,TimeEnum); @@ -2199,11 +2205,11 @@ /*Allocate vector*/ - Vector *vector_control = new Vector(this->vertices->NumberOfVertices()); - Vector *vector_gradient = new Vector(this->vertices->NumberOfVertices()); + Vector *vector_control = new Vector(M[i]*N[i]); + Vector *vector_gradient = new Vector(M[i]*N[i]); /*Fill in vector*/ for(int j=0;jSize();j++){ Element* element=(Element*)elements->GetObjectByOffset(j); - element->ControlToVectors(vector_control,vector_gradient,control_enum); + element->ControlToVectors(vector_control,vector_gradient,control_enum,control_interp[i]); } vector_control->Assemble(); @@ -2216,4 +2222,7 @@ /*Clean up and return*/ xDelete(control_type); + xDelete(control_interp); + xDelete(M); + xDelete(N); } /*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Inputs2/ControlInput2.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs2/ControlInput2.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Inputs2/ControlInput2.cpp (revision 25317) @@ -29,6 +29,4 @@ this->control_id = id; this->layout_enum = input_layout_enum; - - _assert_(interp==P1Enum); switch(this->layout_enum){ Index: /issm/trunk-jpl/src/c/classes/Inputs2/Inputs2.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs2/Inputs2.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Inputs2/Inputs2.cpp (revision 25317) @@ -912,4 +912,5 @@ /*Create it if necessary*/ if(this->inputs[id]){ + int* tenp = xNew(3); if(this->inputs[id]->ObjectEnum()!=PentaInput2Enum) _error_("cannot add Element values to a "<inputs[id]->ObjectEnum())); } Index: /issm/trunk-jpl/src/c/classes/Params/Parameters.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Params/Parameters.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/classes/Params/Parameters.cpp (revision 25317) @@ -386,5 +386,9 @@ /*cast to "passive"*/ + #ifdef _HAVE_AD_ *pscalar=reCast(intermediary); + #else + *pscalar=intermediary; + #endif } /*}}}*/ @@ -401,9 +405,14 @@ /*Make output passive*/ + #ifdef _HAVE_AD_ IssmPDouble* output = xNew(n); for(int i=0;i(vector[i]); + xDelete(vector); + if(pvec) *pvec = output; + #else + if(pvec) *pvec = vector; + #endif /*assign output pointers*/ - if(pvec) *pvec = output; if(pM) *pM = n; }/*}}}*/ Index: /issm/trunk-jpl/src/c/cores/controlm1qn3_core.cpp =================================================================== --- /issm/trunk-jpl/src/c/cores/controlm1qn3_core.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/cores/controlm1qn3_core.cpp (revision 25317) @@ -46,4 +46,10 @@ double *X = NULL; double *G = NULL; + + /*Get control sizes*/ + int* M = NULL; + int* N = NULL; + femmodel->parameters->FindParam(&M,NULL,ControlInputSizeMEnum); + femmodel->parameters->FindParam(&N,NULL,ControlInputSizeNEnum); /*Recover some parameters*/ @@ -73,11 +79,10 @@ long io = 6; /*Channel number for the output*/ - /*Optimization criterions*/ + /*Optimization criterions (need to be cast to long for m1qn3)*/ long niter = long(maxsteps); /*Maximum number of iterations*/ - long nsim = long(maxiter);/*Maximum number of function calls*/ + long nsim = long(maxiter); /*Maximum number of function calls*/ /*Get initial guess*/ GetPassiveVectorFromControlInputsx(&X,&intn,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,"value"); - _assert_(intn==numberofvertices*num_controls); /*Get problem dimension and initialize gradient and initial guess*/ @@ -86,9 +91,11 @@ /*Scale control for M1QN3*/ + int offset = 0; for(int c=0;c(scaling_factors[c]); } + offset += M[c]*N[c]; } @@ -142,11 +149,14 @@ GetPassiveVectorFromControlInputsx(&XL,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,"lowerbound"); GetPassiveVectorFromControlInputsx(&XU,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,"upperbound"); + + offset = 0; for(int c=0;cXU[index]) X[index]=XU[index]; if(X[index](M); + xDelete(N); xDelete(G); xDelete(X); @@ -201,5 +213,9 @@ int num_responses,num_controls,numberofvertices,solution_type; double* scaling_factors = NULL; + int* M = NULL; + int* N = NULL; femmodel->parameters->FindParam(&num_responses,InversionNumCostFunctionsEnum); + femmodel->parameters->FindParam(&N,NULL,ControlInputSizeNEnum); + femmodel->parameters->FindParam(&M,NULL,ControlInputSizeNEnum); femmodel->parameters->FindParam(&num_controls,InversionNumControlParametersEnum); femmodel->parameters->FindParamAndMakePassive(&scaling_factors,NULL,InversionControlScalingFactorsEnum); @@ -212,12 +228,16 @@ GetPassiveVectorFromControlInputsx(&XL,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,"lowerbound"); GetPassiveVectorFromControlInputsx(&XU,NULL,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,"upperbound"); + + int offset = 0; for(int c=0;cXU[index]) X[index]=XU[index]; if(X[index](*n); @@ -277,7 +297,8 @@ /*Constrain Gradient*/ IssmDouble Gnorm = 0.; + offset = 0; for(int c=0;c=XU[index]) G[index]=0.; if(X[index]<=XL[index]) G[index]=0.; @@ -286,4 +307,5 @@ Gnorm += G[index]*G[index]; } + offset += M[c]*N[c]; } Gnorm = sqrt(Gnorm); @@ -296,4 +318,6 @@ /*Clean-up and return*/ *Jlisti = (*Jlisti) +1; + xDelete(M); + xDelete(N); xDelete(XU); xDelete(XL); Index: /issm/trunk-jpl/src/c/modules/ControlInputSetGradientx/ControlInputSetGradientx.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/ControlInputSetGradientx/ControlInputSetGradientx.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/modules/ControlInputSetGradientx/ControlInputSetGradientx.cpp (revision 25317) @@ -9,50 +9,32 @@ void ControlInputSetGradientx(Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,IssmDouble* gradient){ - bool isautodiff; - parameters->FindParam(&isautodiff,AutodiffIsautodiffEnum); - if(isautodiff){ - /*Intermediaries*/ - int num_controls; - int *control_type = NULL; - int* M_all = NULL; - int* N_all = NULL; + /*Intermediaries*/ + int num_controls; + int *control_type = NULL; + int* M_all = NULL; + int* N_all = NULL; + int* interp_all = NULL; - /*Retrieve some parameters*/ - parameters->FindParam(&num_controls,InversionNumControlParametersEnum); - parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); - parameters->FindParam(&M_all,NULL,ControlInputSizeMEnum); - parameters->FindParam(&N_all,NULL,ControlInputSizeNEnum); + /*Retrieve some parameters*/ + parameters->FindParam(&num_controls,InversionNumControlParametersEnum); + parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); + parameters->FindParam(&M_all,NULL,ControlInputSizeMEnum); + parameters->FindParam(&N_all,NULL,ControlInputSizeNEnum); + parameters->FindParam(&interp_all,NULL,ControlInputInterpolationEnum); - int offset = 0; - for(int i=0;iSize();j++){ - Element* element=(Element*)elements->GetObjectByOffset(j); - element->ControlInputSetGradient(gradient,control_type[i],i,offset,N_all[i],M_all[i]); - } - offset+=M_all[i]*N_all[i]; + int offset = 0; + for(int i=0;iSize();j++){ + Element* element=(Element*)elements->GetObjectByOffset(j); + element->ControlInputSetGradient(gradient,control_type[i],i,offset,M_all[i],N_all[i],interp_all[i]); } + offset+=M_all[i]*N_all[i]; + } - /*Clean up and return*/ - xDelete(control_type); - } - else{ - int num_controls; - int *control_type = NULL; - - /*Retrieve some parameters*/ - parameters->FindParam(&num_controls,InversionNumControlParametersEnum); - parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); - - int offset = 0; - for(int i=0;iSize();j++){ - Element* element=(Element*)elements->GetObjectByOffset(j); - element->ControlInputSetGradient(gradient,control_type[i],i); - } - } - - /*Clean up and return*/ - xDelete(control_type); - } + /*Clean up and return*/ + xDelete(control_type); + xDelete(M_all); + xDelete(N_all); + xDelete(interp_all); } @@ -60,6 +42,5 @@ /*Serialize gradient*/ - IssmDouble* serial_gradient=NULL; - serial_gradient=gradient->ToMPISerial(); + IssmDouble* serial_gradient=gradient->ToMPISerial(); ControlInputSetGradientx(elements,nodes,vertices, loads, materials, parameters,serial_gradient); Index: /issm/trunk-jpl/src/c/modules/GetVectorFromControlInputsx/GetVectorFromControlInputsx.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/GetVectorFromControlInputsx/GetVectorFromControlInputsx.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/modules/GetVectorFromControlInputsx/GetVectorFromControlInputsx.cpp (revision 25317) @@ -9,72 +9,38 @@ void GetVectorFromControlInputsx(Vector** pvector, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters,const char* data){/*{{{*/ - bool isautodiff; - parameters->FindParam(&isautodiff,AutodiffIsautodiffEnum); - if(isautodiff){ - int* N = NULL; - int* M = NULL; - int num_controls; - int* control_type = NULL; - Vector* vector=NULL; + int num_controls; + int* N = NULL; + int* M = NULL; + int* control_type = NULL; - /*Retrieve some parameters*/ - parameters->FindParam(&num_controls,InversionNumControlParametersEnum); - parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); - parameters->FindParam(&N,NULL,ControlInputSizeNEnum); - parameters->FindParam(&M,NULL,ControlInputSizeMEnum); + /*Retrieve some parameters*/ + parameters->FindParam(&num_controls,InversionNumControlParametersEnum); + parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); + parameters->FindParam(&N,NULL,ControlInputSizeNEnum); + parameters->FindParam(&M,NULL,ControlInputSizeMEnum); - /*1. Get vector size*/ - int size = 0; - for(int i=0;i(size); + /*2. Allocate vector*/ + Vector* vector=new Vector(size); - /*3. Populate vector*/ - int offset = 0; - for(int i=0;iSize();j++){ - Element* element=(Element*)elements->GetObjectByOffset(j); - element->GetVectorFromControlInputs(vector,control_type[i],i,data,offset); - } - offset += M[i]*N[i]; + /*3. Populate vector*/ + int offset = 0; + for(int i=0;iSize();j++){ + Element* element=(Element*)elements->GetObjectByOffset(j); + element->GetVectorFromControlInputs(vector,control_type[i],i,data,offset); } + offset += M[i]*N[i]; + } + vector->Assemble(); - vector->Assemble(); - - /*Assign output pointers:*/ - xDelete(control_type); - xDelete(M); - xDelete(N); - - *pvector=vector; - } - else{ - int num_controls; - int* control_type = NULL; - Vector* vector=NULL; - - /*Retrieve some parameters*/ - parameters->FindParam(&num_controls,InversionNumControlParametersEnum); - parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); - - /*2. Allocate vector*/ - vector=new Vector(num_controls*vertices->NumberOfVertices()); - - /*3. Populate vector*/ - int offset = 0; - for(int i=0;iSize();j++){ - Element* element=(Element*)elements->GetObjectByOffset(j); - element->GetVectorFromControlInputs(vector,control_type[i],i,data); - } - } - vector->Assemble(); - - /*Assign output pointers:*/ - xDelete(control_type); - *pvector=vector; - } - + /*Assign output pointers:*/ + xDelete(control_type); + xDelete(M); + xDelete(N); + *pvector=vector; }/*}}}*/ void GetVectorFromControlInputsx( IssmDouble** pvector,int *pN, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters, const char* data){/*{{{*/ Index: /issm/trunk-jpl/src/c/modules/Gradjx/Gradjx.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/Gradjx/Gradjx.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/modules/Gradjx/Gradjx.cpp (revision 25317) @@ -9,9 +9,12 @@ void Gradjx(Vector** pgradient,IssmDouble** pnorm_list, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters){ - int numberofvertices; - int num_controls,analysisenum; - IssmDouble norm_inf; - IssmDouble *norm_list = NULL; - int *control_type = NULL; + int numberofvertices; + int num_controls,analysisenum; + IssmDouble norm_inf; + IssmDouble *norm_list = NULL; + int *control_type = NULL; + int *control_interp = NULL; + int *M = NULL; + int *N = NULL; Vector *gradient = NULL; Vector **gradient_list = NULL; @@ -21,6 +24,9 @@ /*retrieve some parameters: */ - parameters->FindParam(&num_controls,InversionNumControlParametersEnum); _assert_(num_controls); + parameters->FindParam(&num_controls,InversionNumControlParametersEnum); _assert_(num_controls); parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); + parameters->FindParam(&control_interp,NULL,ControlInputInterpolationEnum); + parameters->FindParam(&M,NULL,ControlInputSizeMEnum); + parameters->FindParam(&N,NULL,ControlInputSizeNEnum); numberofvertices=vertices->NumberOfVertices(); @@ -31,9 +37,8 @@ /*Allocate gradient_list */ gradient_list = xNew*>(num_controls); - norm_list = xNew(num_controls); - for(int i=0;i(num_controls*numberofvertices); - } - gradient=new Vector(num_controls*numberofvertices); + norm_list = xNew(num_controls); + int totalsize = 0; + for(int i=0;i(totalsize); /*Compute all gradient_list*/ @@ -41,5 +46,5 @@ for(int j=0;jSize();j++){ Element* element=(Element*)elements->GetObjectByOffset(j); - analysis->GradientJ(gradient_list[i],element,control_type[i],i); + analysis->GradientJ(gradient_list[i],element,control_type[i],control_interp[i],i); } gradient_list[i]->Assemble(); @@ -48,4 +53,5 @@ /*Add all gradient_list together*/ + gradient=new Vector(totalsize); for(int i=0;iAXPY(gradient_list[i],1.0); @@ -62,4 +68,7 @@ xDelete*>(gradient_list); xDelete(control_type); + xDelete(control_interp); + xDelete(M); + xDelete(N); if(pnorm_list){ *pnorm_list=norm_list; @@ -73,17 +82,13 @@ void Gradjx(IssmDouble** pgradient,IssmDouble** pnorm_list, Elements* elements,Nodes* nodes, Vertices* vertices, Loads* loads, Materials* materials, Parameters* parameters){ - /*output: */ - IssmDouble* gradient=NULL; + /*Get gradient: */ + Vector* vec_gradient=NULL; + Gradjx(&vec_gradient,pnorm_list,elements,nodes, vertices,loads,materials,parameters); - /*intermediary: */ - Vector* vec_gradient=NULL; + /*Serialize*/ + IssmDouble* gradient=vec_gradient->ToMPISerial(); - Gradjx(&vec_gradient,pnorm_list,elements,nodes, vertices,loads,materials,parameters); - gradient=vec_gradient->ToMPISerial(); - - /*Free ressources:*/ + /*Free ressources and assign output pointer*/ delete vec_gradient; - - /*Assign output pointers:*/ *pgradient=gradient; } Index: /issm/trunk-jpl/src/c/modules/ModelProcessorx/Control/UpdateElementsAndMaterialsControl.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/ModelProcessorx/Control/UpdateElementsAndMaterialsControl.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/modules/ModelProcessorx/Control/UpdateElementsAndMaterialsControl.cpp (revision 25317) @@ -110,4 +110,9 @@ } + int* M_all = xNew(num_controls); + int* N_all = xNew(num_controls); + int* Interp_all = xNew(num_controls); + + int offset = 0; for(int i=0;inumberofvertices && N!=1) _error_("not supported yet"); + + /*Transient independents not supported outside of AD*/ + if(N!=1) _error_("Transient controls not supported yet"); + N_all[i] = N; + M_all[i] = M; + + if(M==iomodel->numberofvertices){ + Interp_all[i] = P1Enum; + } + else if(M==iomodel->numberofelements){ + Interp_all[i] = P0Enum; + } + else{ + _error_("Control size not supported"); + } /*Special case if 3d*/ @@ -143,11 +162,18 @@ for(int j=0;jSize();j++){ Element* element=xDynamicCast(elements->GetObjectByOffset(j)); - element->ControlInputCreate(independent,&independents_min[i*iomodel->numberofvertices],&independents_max[i*iomodel->numberofvertices],inputs2,iomodel,M,N,scale,control,i+1); + element->ControlInputCreate(independent,&independents_min[offset],&independents_max[offset],inputs2,iomodel,M,N,scale,control,i+1); } xDelete(independent); - } + + offset += M*N; + } + parameters->AddObject(new IntVecParam(ControlInputSizeMEnum,M_all,num_controls)); + parameters->AddObject(new IntVecParam(ControlInputSizeNEnum,N_all,num_controls)); + parameters->AddObject(new IntVecParam(ControlInputInterpolationEnum,Interp_all,num_controls)); + xDelete(M_all); + xDelete(N_all); + xDelete(Interp_all); xDelete(independents_min); xDelete(independents_max); - /*Free data: */ Index: /issm/trunk-jpl/src/c/modules/SetControlInputsFromVectorx/SetControlInputsFromVectorx.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/SetControlInputsFromVectorx/SetControlInputsFromVectorx.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/modules/SetControlInputsFromVectorx/SetControlInputsFromVectorx.cpp (revision 25317) @@ -9,46 +9,27 @@ void SetControlInputsFromVectorx(FemModel* femmodel,IssmDouble* vector){ - bool isautodiff; - femmodel->parameters->FindParam(&isautodiff,AutodiffIsautodiffEnum); - if(isautodiff){ - int num_controls; - int* control_type = NULL; - int* M = NULL; - int* N = NULL; + int num_controls; + int* control_type = NULL; + int* M = NULL; + int* N = NULL; - /*Retrieve some parameters*/ - femmodel->parameters->FindParam(&num_controls,InversionNumControlParametersEnum); - femmodel->parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); - femmodel->parameters->FindParam(&M,NULL,ControlInputSizeMEnum); - femmodel->parameters->FindParam(&N,NULL,ControlInputSizeNEnum); + /*Retrieve some parameters*/ + femmodel->parameters->FindParam(&num_controls,InversionNumControlParametersEnum); + femmodel->parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); + femmodel->parameters->FindParam(&M,NULL,ControlInputSizeMEnum); + femmodel->parameters->FindParam(&N,NULL,ControlInputSizeNEnum); - int offset = 0; - for(int i=0;ielements->Size();j++){ - Element* element=(Element*)femmodel->elements->GetObjectByOffset(j); - element->SetControlInputsFromVector(vector,control_type[i],i,offset,N[i],M[i]); - } + int offset = 0; + for(int i=0;ielements->Size();j++){ + Element* element=(Element*)femmodel->elements->GetObjectByOffset(j); + element->SetControlInputsFromVector(vector,control_type[i],i,offset,M[i],N[i]); + } offset += M[i]*N[i]; - } + } - xDelete(control_type); - xDelete(M); - xDelete(N); - } - else{ - - int num_controls; - int* control_type = NULL; - femmodel->parameters->FindParam(&num_controls,InversionNumControlParametersEnum); - femmodel->parameters->FindParam(&control_type,NULL,InversionControlParametersEnum); - int offset = 0; - for(int i=0;ielements->Size();j++){ - Element* element=(Element*)femmodel->elements->GetObjectByOffset(j); - element->SetControlInputsFromVector(vector,control_type[i],i); - } - } - xDelete(control_type); - } + xDelete(control_type); + xDelete(M); + xDelete(N); } Index: /issm/trunk-jpl/src/c/shared/Enum/Enum.vim =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/Enum.vim (revision 25316) +++ /issm/trunk-jpl/src/c/shared/Enum/Enum.vim (revision 25317) @@ -111,4 +111,5 @@ syn keyword cConstant ControlInputSizeMEnum syn keyword cConstant ControlInputSizeNEnum +syn keyword cConstant ControlInputInterpolationEnum syn keyword cConstant DamageC1Enum syn keyword cConstant DamageC2Enum Index: /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 25316) +++ /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 25317) @@ -105,4 +105,5 @@ ControlInputSizeMEnum, ControlInputSizeNEnum, + ControlInputInterpolationEnum, DamageC1Enum, DamageC2Enum, Index: /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 25317) @@ -113,4 +113,5 @@ case ControlInputSizeMEnum : return "ControlInputSizeM"; case ControlInputSizeNEnum : return "ControlInputSizeN"; + case ControlInputInterpolationEnum : return "ControlInputInterpolation"; case DamageC1Enum : return "DamageC1"; case DamageC2Enum : return "DamageC2"; Index: /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 25316) +++ /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 25317) @@ -113,4 +113,5 @@ else if (strcmp(name,"ControlInputSizeM")==0) return ControlInputSizeMEnum; else if (strcmp(name,"ControlInputSizeN")==0) return ControlInputSizeNEnum; + else if (strcmp(name,"ControlInputInterpolation")==0) return ControlInputInterpolationEnum; else if (strcmp(name,"DamageC1")==0) return DamageC1Enum; else if (strcmp(name,"DamageC2")==0) return DamageC2Enum; @@ -136,9 +137,9 @@ else if (strcmp(name,"DslComputeFingerprints")==0) return DslComputeFingerprintsEnum; else if (strcmp(name,"EarthId")==0) return EarthIdEnum; - else if (strcmp(name,"EplZigZagCounter")==0) return EplZigZagCounterEnum; else stage=2; } if(stage==2){ - if (strcmp(name,"EsaHElastic")==0) return EsaHElasticEnum; + if (strcmp(name,"EplZigZagCounter")==0) return EplZigZagCounterEnum; + else if (strcmp(name,"EsaHElastic")==0) return EsaHElasticEnum; else if (strcmp(name,"EsaHemisphere")==0) return EsaHemisphereEnum; else if (strcmp(name,"EsaRequestedOutputs")==0) return EsaRequestedOutputsEnum; @@ -259,9 +260,9 @@ else if (strcmp(name,"MasstransportRequestedOutputs")==0) return MasstransportRequestedOutputsEnum; else if (strcmp(name,"MasstransportStabilization")==0) return MasstransportStabilizationEnum; - else if (strcmp(name,"MaterialsBeta")==0) return MaterialsBetaEnum; else stage=3; } if(stage==3){ - if (strcmp(name,"MaterialsEarthDensity")==0) return MaterialsEarthDensityEnum; + if (strcmp(name,"MaterialsBeta")==0) return MaterialsBetaEnum; + else if (strcmp(name,"MaterialsEarthDensity")==0) return MaterialsEarthDensityEnum; else if (strcmp(name,"MaterialsEffectiveconductivityAveraging")==0) return MaterialsEffectiveconductivityAveragingEnum; else if (strcmp(name,"MaterialsHeatcapacity")==0) return MaterialsHeatcapacityEnum; @@ -382,9 +383,9 @@ else if (strcmp(name,"SmbIsfirnwarming")==0) return SmbIsfirnwarmingEnum; else if (strcmp(name,"SmbIsgraingrowth")==0) return SmbIsgraingrowthEnum; - else if (strcmp(name,"SmbIsmelt")==0) return SmbIsmeltEnum; else stage=4; } if(stage==4){ - if (strcmp(name,"SmbIsmungsm")==0) return SmbIsmungsmEnum; + if (strcmp(name,"SmbIsmelt")==0) return SmbIsmeltEnum; + else if (strcmp(name,"SmbIsmungsm")==0) return SmbIsmungsmEnum; else if (strcmp(name,"SmbIsprecipscaled")==0) return SmbIsprecipscaledEnum; else if (strcmp(name,"SmbIssetpddfac")==0) return SmbIssetpddfacEnum; @@ -505,9 +506,9 @@ else if (strcmp(name,"BasalforcingsPicoSubShelfOceanSalinity")==0) return BasalforcingsPicoSubShelfOceanSalinityEnum; else if (strcmp(name,"BasalforcingsPicoSubShelfOceanTemp")==0) return BasalforcingsPicoSubShelfOceanTempEnum; - else if (strcmp(name,"BasalStressx")==0) return BasalStressxEnum; else stage=5; } if(stage==5){ - if (strcmp(name,"BasalStressy")==0) return BasalStressyEnum; + if (strcmp(name,"BasalStressx")==0) return BasalStressxEnum; + else if (strcmp(name,"BasalStressy")==0) return BasalStressyEnum; else if (strcmp(name,"BasalStress")==0) return BasalStressEnum; else if (strcmp(name,"Base")==0) return BaseEnum; @@ -628,9 +629,9 @@ else if (strcmp(name,"HydrologydcSedimentTransmitivity")==0) return HydrologydcSedimentTransmitivityEnum; else if (strcmp(name,"HydrologyDrainageRate")==0) return HydrologyDrainageRateEnum; - else if (strcmp(name,"HydrologyEnglacialInput")==0) return HydrologyEnglacialInputEnum; else stage=6; } if(stage==6){ - if (strcmp(name,"HydrologyGapHeight")==0) return HydrologyGapHeightEnum; + if (strcmp(name,"HydrologyEnglacialInput")==0) return HydrologyEnglacialInputEnum; + else if (strcmp(name,"HydrologyGapHeight")==0) return HydrologyGapHeightEnum; else if (strcmp(name,"HydrologyHead")==0) return HydrologyHeadEnum; else if (strcmp(name,"HydrologyHeadOld")==0) return HydrologyHeadOldEnum; @@ -751,9 +752,9 @@ else if (strcmp(name,"SmbDzTop")==0) return SmbDzTopEnum; else if (strcmp(name,"SmbDzini")==0) return SmbDziniEnum; - else if (strcmp(name,"SmbEAir")==0) return SmbEAirEnum; else stage=7; } if(stage==7){ - if (strcmp(name,"SmbEC")==0) return SmbECEnum; + if (strcmp(name,"SmbEAir")==0) return SmbEAirEnum; + else if (strcmp(name,"SmbEC")==0) return SmbECEnum; else if (strcmp(name,"SmbECDt")==0) return SmbECDtEnum; else if (strcmp(name,"SmbECini")==0) return SmbECiniEnum; @@ -874,9 +875,9 @@ else if (strcmp(name,"VyMesh")==0) return VyMeshEnum; else if (strcmp(name,"VyObs")==0) return VyObsEnum; - else if (strcmp(name,"Vz")==0) return VzEnum; else stage=8; } if(stage==8){ - if (strcmp(name,"VzFS")==0) return VzFSEnum; + if (strcmp(name,"Vz")==0) return VzEnum; + else if (strcmp(name,"VzFS")==0) return VzFSEnum; else if (strcmp(name,"VzHO")==0) return VzHOEnum; else if (strcmp(name,"VzMesh")==0) return VzMeshEnum; @@ -997,9 +998,9 @@ else if (strcmp(name,"AggressiveMigration")==0) return AggressiveMigrationEnum; else if (strcmp(name,"AmrBamg")==0) return AmrBamgEnum; - else if (strcmp(name,"AmrNeopz")==0) return AmrNeopzEnum; else stage=9; } if(stage==9){ - if (strcmp(name,"AndroidFrictionCoefficient")==0) return AndroidFrictionCoefficientEnum; + if (strcmp(name,"AmrNeopz")==0) return AmrNeopzEnum; + else if (strcmp(name,"AndroidFrictionCoefficient")==0) return AndroidFrictionCoefficientEnum; else if (strcmp(name,"Arrhenius")==0) return ArrheniusEnum; else if (strcmp(name,"AutodiffJacobian")==0) return AutodiffJacobianEnum; @@ -1120,9 +1121,9 @@ else if (strcmp(name,"Gset")==0) return GsetEnum; else if (strcmp(name,"Gsl")==0) return GslEnum; - else if (strcmp(name,"HOApproximation")==0) return HOApproximationEnum; else stage=10; } if(stage==10){ - if (strcmp(name,"HOFSApproximation")==0) return HOFSApproximationEnum; + if (strcmp(name,"HOApproximation")==0) return HOApproximationEnum; + else if (strcmp(name,"HOFSApproximation")==0) return HOFSApproximationEnum; else if (strcmp(name,"Hook")==0) return HookEnum; else if (strcmp(name,"HydrologyDCEfficientAnalysis")==0) return HydrologyDCEfficientAnalysisEnum; @@ -1243,9 +1244,9 @@ else if (strcmp(name,"P0Array")==0) return P0ArrayEnum; else if (strcmp(name,"P0DG")==0) return P0DGEnum; - else if (strcmp(name,"P1DG")==0) return P1DGEnum; else stage=11; } if(stage==11){ - if (strcmp(name,"P1P1")==0) return P1P1Enum; + if (strcmp(name,"P1DG")==0) return P1DGEnum; + else if (strcmp(name,"P1P1")==0) return P1P1Enum; else if (strcmp(name,"P1P1GLS")==0) return P1P1GLSEnum; else if (strcmp(name,"P1bubble")==0) return P1bubbleEnum; @@ -1366,9 +1367,9 @@ else if (strcmp(name,"BalancethicknessDiffusionCoefficient")==0) return BalancethicknessDiffusionCoefficientEnum; else if (strcmp(name,"BilinearInterp")==0) return BilinearInterpEnum; - else if (strcmp(name,"CalvingdevCoeff")==0) return CalvingdevCoeffEnum; else stage=12; } if(stage==12){ - if (strcmp(name,"DeviatoricStress")==0) return DeviatoricStressEnum; + if (strcmp(name,"CalvingdevCoeff")==0) return CalvingdevCoeffEnum; + else if (strcmp(name,"DeviatoricStress")==0) return DeviatoricStressEnum; else if (strcmp(name,"EtaAbsGradient")==0) return EtaAbsGradientEnum; else if (strcmp(name,"MeshZ")==0) return MeshZEnum;