Changeset 4047


Ignore:
Timestamp:
06/15/10 11:13:28 (15 years ago)
Author:
Eric.Larour
Message:

Starting cleanup of control_core.
Introduced gradient_core and adjoint_core solutions. Trying to get the gradient and the adjoint
inside the elements, instead of as node vectors.

Location:
issm/trunk/src
Files:
7 added
2 deleted
86 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk/src/c/DataSet/DataSet.cpp

    r4043 r4047  
    893893                        node=(Node*)(*object);
    894894
    895                         if (node->InAnalysis(analysis_type)){
     895                        /*Check that this node corresponds to our analysis currently being carried out: */
    896896
    897897                                /*Ok, this object is a node, ask it to plug values into partition: */

  • issm/trunk/src/c/EnumDefinitions/EnumAsString.cpp

    r4043 r4047  
    173173                case PressureOldEnum : return "PressureOld";
    174174                case QmuPressureEnum : return "QmuPressure";
    175                 case StokesPressureEnum : return "StokesPressure";
     175                case PressureStokesEnum : return "PressureStokes";
    176176                case ResetPenaltiesEnum : return "ResetPenalties";
    177177                case RheologyBEnum : return "RheologyB";
    … …  
    217217                case CmMaxEnum : return "CmMax";
    218218                case CmMinEnum : return "CmMin";
     219                case GradientEnum : return "Gradient";
    219220                case ConnectivityEnum : return "Connectivity";
    220221                case ControlParameterEnum : return "ControlParameter";

  • issm/trunk/src/c/EnumDefinitions/EnumDefinitions.h

    r4046 r4047  
    252252        CmMaxEnum,
    253253        CmMinEnum,
     254        GradientEnum,
    254255        ConnectivityEnum,
    255256        ControlParameterEnum,

  • issm/trunk/src/c/EnumDefinitions/StringAsEnum.cpp

    r4043 r4047  
    171171        else if (strcmp(name,"PressureOld")==0) return PressureOldEnum;
    172172        else if (strcmp(name,"QmuPressure")==0) return QmuPressureEnum;
    173         else if (strcmp(name,"StokesPressure")==0) return StokesPressureEnum;
     173        else if (strcmp(name,"PressureStokes")==0) return PressureStokesEnum;
    174174        else if (strcmp(name,"ResetPenalties")==0) return ResetPenaltiesEnum;
    175175        else if (strcmp(name,"RheologyB")==0) return RheologyBEnum;
    … …  
    215215        else if (strcmp(name,"CmMax")==0) return CmMaxEnum;
    216216        else if (strcmp(name,"CmMin")==0) return CmMinEnum;
     217        else if (strcmp(name,"Gradient")==0) return GradientEnum;
    217218        else if (strcmp(name,"Connectivity")==0) return ConnectivityEnum;
    218219        else if (strcmp(name,"ControlParameter")==0) return ControlParameterEnum;

  • issm/trunk/src/c/Makefile.am

    r4046 r4047  
    391391                                        ./modules/DuplicateInputx/DuplicateInputx.h\
    392392                                        ./modules/DuplicateInputx/DuplicateInputx.cpp\
     393                                        ./modules/ScaleInputx/ScaleInputx.h\
     394                                        ./modules/ScaleInputx/ScaleInputx.cpp\
    393395                                        ./modules/ControlConstrainx/ControlConstrainx.h\
    394396                                        ./modules/ControlConstrainx/ControlConstrainx.cpp\
    … …  
    887889                                        ./modules/DuplicateInputx/DuplicateInputx.h\
    888890                                        ./modules/DuplicateInputx/DuplicateInputx.cpp\
     891                                        ./modules/ScaleInputx/ScaleInputx.h\
     892                                        ./modules/ScaleInputx/ScaleInputx.cpp\
    889893                                        ./modules/ControlConstrainx/ControlConstrainx.h\
    890894                                        ./modules/ControlConstrainx/ControlConstrainx.cpp\

  • issm/trunk/src/c/modules/Gradjx/Gradjx.cpp

    r4043 r4047  
    1010#include "../../EnumDefinitions/EnumDefinitions.h"
    1111
    12 void Gradjx( Vec* pgrad_g, int numberofnodes, DataSet* elements,DataSet* nodes, DataSet* vertices, DataSet* loads, DataSet* materials, Parameters* parameters,int control_type){
     12void Gradjx( DataSet* elements,DataSet* nodes, DataSet* vertices, DataSet* loads, DataSet* materials, Parameters* parameters,int control_type){
    1313
    14         /*Intermediary*/
    1514        int i;
    16         Element* element=NULL;
     15        int dim;
    1716
    18         /*output: */
    19         Vec grad_g=NULL;
    2017
    2118        /*First, get elements and loads configured: */
    … …  
    2421        parameters->Configure(elements,loads, nodes,vertices, materials,parameters);
    2522
    26         /*Allocate grad_g: */
    27         grad_g=NewVec(numberofnodes);
     23        /*retrieve some parameters: */
     24        parameters->FindParam(&dim,DimEnum);
     25
    2826
    2927        /*Compute gradients: */
    3028        for (i=0;i<elements->Size();i++){
    31                 element=(Element*)elements->GetObjectByOffset(i);
    32                 element->Gradj(grad_g,control_type);
     29                Element* element=(Element*)elements->GetObjectByOffset(i);
     30                element->Gradj(control_type);
    3331        }
    3432
    35         /*Assemble vector: */
    36         VecAssemblyBegin(grad_g);
    37         VecAssemblyEnd(grad_g);
     33        /*Extrude if needed: */
     34        if(dim==3) InputsExtrudex( elements,nodes,vertices,loads,materials,parameters,GradientEnum,0);
    3835
    39         /*Assign output pointers: */
    40         *pgrad_g=grad_g;
    4136}

  • issm/trunk/src/c/modules/Gradjx/Gradjx.h

    r3913 r4047  
    1010
    1111/* local prototypes: */
    12 void Gradjx( Vec* pgrad_g, int numberofnodes, DataSet* elements,DataSet* nodes, DataSet* vertices, DataSet* loads, DataSet* materials,  Parameters* parameters,
    13                         int analysis_type,int sub_analysis_type,int control_type);
     12void Gradjx( DataSet* elements,DataSet* nodes, DataSet* vertices, DataSet* loads, DataSet* materials,  Parameters* parameters, int control_type);
    1413
    1514#endif  /* _GRADJX_H */

  • issm/trunk/src/c/modules/modules.h

    r4046 r4047  
    7979#include "./MaxAbsVzx/MaxAbsVzx.h"
    8080#include "./DuplicateInputx/DuplicateInputx.h"
     81#include "./ScaleInputx/ScaleInputx.h"
    8182
    8283#endif

  • issm/trunk/src/c/objects/Elements/Beam.cpp

    r4046 r4047  
    562562/*}}}*/
    563563/*FUNCTION Beam::Gradj{{{1*/
    564 void  Beam::Gradj(Vec, int control_type){
     564void  Beam::Gradj(int control_type){
    565565        ISSMERROR(" not supported yet!");
    566566}
    567567/*}}}*/
    568568/*FUNCTION Beam::GradjB{{{1*/
    569 void  Beam::GradjB(Vec){
     569void  Beam::GradjB(void){
    570570        ISSMERROR(" not supported yet!");
    571571}
    572572/*}}}*/
    573573/*FUNCTION Beam::GradjDrag{{{1*/
    574 void  Beam::GradjDrag(Vec){
     574void  Beam::GradjDrag(void){
    575575        ISSMERROR(" not supported yet!");
    576576}
    … …  
    990990}
    991991/*}}}*/
     992/*FUNCTION Beam::ScaleInput(int enum_type,double scale_factor){{{1*/
     993void  Beam::ScaleInput(int enum_type,double scale_factor){
     994
     995        Input* input=NULL;
     996
     997        /*Make a copy of the original input: */
     998        input=(Input*)this->inputs->GetInput(enum_type);
     999
     1000        /*Scale: */
     1001        input->Scale(scale_factor);
     1002}
     1003/*}}}*/

  • issm/trunk/src/c/objects/Elements/Beam.h

    r4046 r4047  
    9595                void  MaxAbsVz(double* pmaxabsvz, bool process_units);
    9696                void  DuplicateInput(int original_enum,int new_enum);
     97                void  ScaleInput(int enum_type,double scale_factor);
    9798
    9899                /*}}}*/
    … …  
    103104                void  GetThicknessList(double* thickness_list);
    104105                void  Du(Vec);
    105                 void  Gradj(Vec,  int control_type);
    106                 void  GradjDrag(Vec);
    107                 void  GradjB(Vec);
     106                void  Gradj(int control_type);
     107                void  GradjDrag(void);
     108                void  GradjB(void);
    108109                double Misfit(void);
    109110                double SurfaceArea(void);

  • issm/trunk/src/c/objects/Elements/Element.h

    r4046 r4047  
    3737                virtual void   GetBedList(double* bed_list)=0;
    3838                virtual void   Du(Vec du_g)=0;
    39                 virtual void   Gradj(Vec grad_g,int control_type)=0;
    40                 virtual void   GradjDrag(Vec grad_g)=0;
    41                 virtual void   GradjB(Vec grad_g)=0;
     39                virtual void   Gradj(int control_type)=0;
     40                virtual void   GradjDrag(void)=0;
     41                virtual void   GradjB(void)=0;
    4242                virtual double Misfit(void)=0;
    4343                virtual double CostFunction(void)=0;
    … …  
    6565                virtual void   MaxAbsVz(double* pmaxabsvz, bool process_units)=0;
    6666                virtual void   DuplicateInput(int original_enum,int new_enum)=0;
     67                virtual void   ScaleInput(int enum_type,double scale_factor)=0;
    6768
    6869

  • issm/trunk/src/c/objects/Elements/Penta.cpp

    r4046 r4047  
    45554555/*}}}*/
    45564556/*FUNCTION Penta::Gradj {{{1*/
    4557 void  Penta::Gradj(Vec grad_g,int control_type){
     4557void  Penta::Gradj(int control_type){
    45584558
    45594559        /*inputs: */
    … …  
    45674567
    45684568        if (control_type==DragCoefficientEnum){
    4569                 GradjDrag( grad_g);
     4569                GradjDrag();
    45704570        }
    45714571        else if (control_type=RheologyBEnum){
    4572                 GradjB( grad_g);
     4572                GradjB();
    45734573        }
    45744574        else ISSMERROR("%s%i","control type not supported yet: ",control_type);
    … …  
    45764576/*}}}*/
    45774577/*FUNCTION Penta::GradjDrag {{{1*/
    4578 void  Penta::GradjDrag(Vec grad_g){
     4578void  Penta::GradjDrag(void){
    45794579
    45804580        Tria* tria=NULL;
     4581        TriaVertexInput* triavertexinput=NULL;
     4582        double gradient[6]={0,0,0,0,0,0};
    45814583
    45824584        /*inputs: */
    … …  
    45864588        int analysis_type,sub_analysis_type;
    45874589
    4588         /*retrive parameters: */
     4590        /*retrieve parameters: */
    45894591        parameters->FindParam(&analysis_type,AnalysisTypeEnum);
    45904592        parameters->FindParam(&sub_analysis_type,AnalysisTypeEnum);
    … …  
    46084610                /*MacAyeal or Pattyn*/
    46094611                tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria.
    4610                 tria->GradjDrag( grad_g);
    4611                 delete tria;
    4612                 return;
     4612                tria->GradjDrag();
    46134613        }
    46144614        else if (sub_analysis_type==StokesAnalysisEnum){
    … …  
    46164616                /*Stokes*/
    46174617                tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria.
    4618                 tria->GradjDragStokes( grad_g);
    4619                 delete tria;
    4620                 return;
     4618                tria->GradjDragStokes();
    46214619        }
    46224620        else ISSMERROR("%s%i%s\n","sub_analysis: ",sub_analysis_type," not supported yet");
     4621
     4622        /*Now, the tria  has a GradientEnum input. Take it, and make a PentaVertexInput out of it, then delete the tria: */
     4623        triavertexinput=tria->inputs->GetInput(GradientEnum);
     4624        for(i=0;i<3;i++)gradient[i]=triavertexinput->values[i];
     4625        this->inputs->AddInput(new PentaVertexInput(GradientEnum,&gradient[0]));
     4626       
     4627        delete tria;
     4628
    46234629}
    46244630/*}}}*/
    46254631/*FUNCTION Penta::GradjB {{{1*/
    4626 void  Penta::GradjB(Vec grad_g){
     4632void  Penta::GradjB(void){
    46274633
    46284634        Tria* tria=NULL;
     4635        TriaVertexInput* triavertexinput=NULL;
     4636        double gradient[6]={0,0,0,0,0,0};
    46294637
    46304638        /*inputs: */
    … …  
    46484656                 * and compute gardj*/
    46494657                tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria (lower face).
    4650                 tria->GradjB(grad_g);
    4651                 delete tria;
    4652                 return;
     4658                tria->GradjB();
    46534659        }
    46544660        else{
    46554661                /*B is a 2d field, use MacAyeal(2d) gradient even if it is Stokes or Pattyn*/
    46564662                tria=(Tria*)SpawnTria(0,1,2); //grids 0, 1 and 2 make the new tria (lower face).
    4657                 tria->GradjB(grad_g);
    4658                 delete tria;
    4659                 return;
    4660         }
     4663                tria->GradjB();
     4664        }
     4665       
     4666        /*Now, the tria  has a GradientEnum input. Take it, and make a PentaVertexInput out of it, then delete the tria: */
     4667        triavertexinput=tria->inputs->GetInput(GradientEnum);
     4668        for(i=0;i<3;i++)gradient[i]=triavertexinput->values[i];
     4669        this->inputs->AddInput(new PentaVertexInput(GradientEnum,&gradient[0]));
     4670       
     4671        delete tria;
     4672
    46614673}
    46624674/*}}}*/
    … …  
    53285340}
    53295341/*}}}*/
     5342/*FUNCTION Penta::ScaleInput(int enum_type,double scale_factor){{{1*/
     5343void  Penta::ScaleInput(int enum_type,double scale_factor){
     5344
     5345        Input* input=NULL;
     5346
     5347        /*Make a copy of the original input: */
     5348        input=(Input*)this->inputs->GetInput(enum_type);
     5349
     5350        /*Scale: */
     5351        input->Scale(scale_factor);
     5352}
     5353/*}}}*/

  • issm/trunk/src/c/objects/Elements/Penta.h

    r4046 r4047  
    8686                bool   GetOnBed();
    8787                void  Du(Vec du_gg);
    88                 void  Gradj(Vec grad_gg,int control_type);
    89                 void  GradjDrag(Vec grad_gg);
    90                 void  GradjB(Vec grad_gg);
     88                void  Gradj(int control_type);
     89                void  GradjDrag(void);
     90                void  GradjB(void);
    9191                double Misfit(void);
    9292                double SurfaceArea(void);
    … …  
    162162                void  MaxAbsVz(double* pmaxabsvz, bool process_units);
    163163                void  DuplicateInput(int original_enum,int new_enum);
     164                void  ScaleInput(int enum_type,double scale_factor);
    164165
    165166

  • issm/trunk/src/c/objects/Elements/Sing.cpp

    r4046 r4047  
    383383/*}}}*/
    384384/*FUNCTION Sing::Gradj {{{1*/
    385 void  Sing::Gradj(Vec,  int control_type){
     385void  Sing::Gradj(int control_type){
    386386        ISSMERROR(" not supported yet!");
    387387}
    388388/*}}}*/
    389389/*FUNCTION Sing::GradB {{{1*/
    390 void  Sing::GradjB(Vec){
     390void  Sing::GradjB(void){
    391391        ISSMERROR(" not supported yet!");
    392392}
    393393/*}}}*/
    394394/*FUNCTION Sing::GradjDrag {{{1*/
    395 void  Sing::GradjDrag(Vec){
     395void  Sing::GradjDrag(void){
    396396        ISSMERROR(" not supported yet!");
    397397}
    … …  
    696696}
    697697/*}}}*/
     698/*FUNCTION Sing::ScaleInput(int enum_type,double scale_factor){{{1*/
     699void  Sing::ScaleInput(int enum_type,double scale_factor){
     700
     701        Input* input=NULL;
     702
     703        /*Make a copy of the original input: */
     704        input=(Input*)this->inputs->GetInput(enum_type);
     705
     706        /*Scale: */
     707        input->Scale(scale_factor);
     708}
     709/*}}}*/

  • issm/trunk/src/c/objects/Elements/Sing.h

    r4046 r4047  
    9494                void  MaxAbsVz(double* pmaxabsvz, bool process_units);
    9595                void  DuplicateInput(int original_enum,int new_enum);
     96                void  ScaleInput(int enum_type,double scale_factor);
    9697
    9798
    … …  
    103104                void  GetThicknessList(double* thickness_list);
    104105                void  Du(Vec);
    105                 void  Gradj(Vec, int control_type);
    106                 void  GradjDrag(Vec);
    107                 void  GradjB(Vec);
     106                void  Gradj(int control_type);
     107                void  GradjDrag(void);
     108                void  GradjB(void);
    108109                double Misfit(void);
    109110                double SurfaceArea(void);

  • issm/trunk/src/c/objects/Elements/Tria.cpp

    r4046 r4047  
    41254125/*}}}*/
    41264126/*FUNCTION Tria::Gradj {{{1*/
    4127 void  Tria::Gradj(Vec grad_g,int control_type){
     4127void  Tria::Gradj(int control_type){
    41284128
    41294129        /*inputs: */
    … …  
    41374137
    41384138        if (control_type==DragCoefficientEnum){
    4139                 GradjDrag( grad_g);
     4139                GradjDrag();
    41404140        }
    41414141        else if (control_type==RheologyBEnum){
    4142                 GradjB( grad_g);
     4142                GradjB();
    41434143        }
    41444144        else ISSMERROR("%s%i","control type not supported yet: ",control_type);
     4145}
     4146/*}}}*/
     4147/*FUNCTION Tria::GradjDrag {{{1*/
     4148void  Tria::GradjDrag(void){
     4149
     4150
     4151        int i;
     4152
     4153        /* node data: */
     4154        const int    numgrids=3;
     4155        const int    NDOF2=2;
     4156        const int    numdof=NDOF2*numgrids;
     4157        double       xyz_list[numgrids][3];
     4158        int          doflist1[numgrids];
     4159        double       dh1dh3[NDOF2][numgrids];
     4160
     4161        /* grid data: */
     4162        double adjx_list[numgrids];
     4163        double adjy_list[numgrids];
     4164
     4165        /* gaussian points: */
     4166        int     num_gauss,ig;
     4167        double* first_gauss_area_coord  =  NULL;
     4168        double* second_gauss_area_coord =  NULL;
     4169        double* third_gauss_area_coord  =  NULL;
     4170        double* gauss_weights           =  NULL;
     4171        double  gauss_weight;
     4172        double  gauss_l1l2l3[3];
     4173
     4174        /* parameters: */
     4175        double  dk[NDOF2];
     4176        double  vx,vy;
     4177        double  lambda,mu;
     4178        double  bed,thickness,Neff;
     4179        double  alpha_complement;
     4180        int     drag_type;
     4181        double  drag;
     4182        Friction* friction=NULL;
     4183
     4184        /*element vector at the gaussian points: */
     4185        double  grade_g[numgrids]={0.0};
     4186        double  grade_g_gaussian[numgrids];
     4187
     4188        /* Jacobian: */
     4189        double Jdet;
     4190
     4191        /*nodal functions: */
     4192        double l1l2l3[3];
     4193
     4194        /* strain rate: */
     4195        double epsilon[3]; /* epsilon=[exx,eyy,exy];*/
     4196
     4197        /*inputs: */
     4198        bool shelf;
     4199
     4200        /*parameters: */
     4201        double  cm_noisedmp;
     4202        double  cm_mindmp_slope;
     4203        double  cm_mindmp_value;
     4204        double  cm_maxdmp_value;
     4205        double  cm_maxdmp_slope;
     4206
     4207        int analysis_type;
     4208
     4209        /*retrive parameters: */
     4210        parameters->FindParam(&analysis_type,AnalysisTypeEnum);
     4211
     4212        /*retrieve inputs :*/
     4213        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
     4214
     4215        /*retrieve some parameters: */
     4216        this->parameters->FindParam(&cm_noisedmp,CmNoiseDmpEnum);
     4217        this->parameters->FindParam(&cm_mindmp_value,CmMinDmpValueEnum);
     4218        this->parameters->FindParam(&cm_mindmp_slope,CmMinDmpSlopeEnum);
     4219        this->parameters->FindParam(&cm_maxdmp_value,CmMaxDmpValueEnum);
     4220        this->parameters->FindParam(&cm_maxdmp_slope,CmMaxDmpSlopeEnum);
     4221
     4222
     4223        /*Get out if shelf*/
     4224        if(shelf)return;
     4225
     4226        /* Get node coordinates and dof list: */
     4227        GetVerticesCoordinates(&xyz_list[0][0], nodes, numgrids);
     4228        GetDofList1(&doflist1[0]);
     4229
     4230        /*Build frictoin element, needed later: */
     4231        inputs->GetParameterValue(&drag_type,DragTypeEnum);
     4232        friction=new Friction("2d",inputs,matpar,analysis_type);
     4233
     4234        /* Get gaussian points and weights (make this a statically initialized list of points? fstd): */
     4235        GaussTria( &num_gauss, &first_gauss_area_coord, &second_gauss_area_coord, &third_gauss_area_coord, &gauss_weights, 4);
     4236
     4237        /* Start  looping on the number of gaussian points: */
     4238        for (ig=0; ig<num_gauss; ig++){
     4239                /*Pick up the gaussian point: */
     4240                gauss_weight=*(gauss_weights+ig);
     4241                gauss_l1l2l3[0]=*(first_gauss_area_coord+ig);
     4242                gauss_l1l2l3[1]=*(second_gauss_area_coord+ig);
     4243                gauss_l1l2l3[2]=*(third_gauss_area_coord+ig);
     4244
     4245                /*Build alpha_complement_list: */
     4246                if (drag_type==2) friction->GetAlphaComplement(&alpha_complement, gauss_l1l2l3,VxAverageEnum,VyAverageEnum);
     4247                else alpha_complement=0;
     4248       
     4249                /*Recover alpha_complement and k: */
     4250                inputs->GetParameterValue(&drag, gauss_l1l2l3,DragCoefficientEnum);
     4251
     4252                /*recover lambda and mu: */
     4253                inputs->GetParameterValue(&lambda, gauss_l1l2l3,AdjointxEnum);
     4254                inputs->GetParameterValue(&mu, gauss_l1l2l3,AdjointyEnum);
     4255                       
     4256                /*recover vx and vy: */
     4257                inputs->GetParameterValue(&vx, gauss_l1l2l3,VxEnum);
     4258                inputs->GetParameterValue(&vy, gauss_l1l2l3,VyEnum);
     4259
     4260                /* Get Jacobian determinant: */
     4261                GetJacobianDeterminant2d(&Jdet, &xyz_list[0][0],gauss_l1l2l3);
     4262               
     4263                /* Get nodal functions value at gaussian point:*/
     4264                GetNodalFunctions(l1l2l3, gauss_l1l2l3);
     4265
     4266                /*Get nodal functions derivatives*/
     4267                GetNodalFunctionsDerivatives(&dh1dh3[0][0],&xyz_list[0][0],gauss_l1l2l3);
     4268
     4269                /*Get k derivative: dk/dx */
     4270                inputs->GetParameterDerivativeValue(&dk[0],&xyz_list[0][0],&gauss_l1l2l3[0],DragCoefficientEnum);
     4271
     4272                /*Build gradje_g_gaussian vector (actually -dJ/ddrag): */
     4273                for (i=0;i<numgrids;i++){
     4274
     4275                        //standard term dJ/dki
     4276                        grade_g_gaussian[i]=-2*drag*alpha_complement*((lambda*vx+mu*vy))*Jdet*gauss_weight*l1l2l3[i];
     4277
     4278                        //noise dampening d/dki(1/2*(dk/dx)^2)
     4279                        grade_g_gaussian[i]+=-cm_noisedmp*Jdet*gauss_weight*(dh1dh3[0][i]*dk[0]+dh1dh3[1][i]*dk[1]);
     4280                       
     4281                        //min dampening
     4282                        if(drag<cm_mindmp_value){
     4283                                grade_g_gaussian[i]+=cm_mindmp_slope*Jdet*gauss_weight*l1l2l3[i];
     4284                        }
     4285
     4286                        //max dampening
     4287                        if(drag>cm_maxdmp_value){
     4288                                grade_g_gaussian[i]+= - cm_maxdmp_slope*Jdet*gauss_weight*l1l2l3[i];
     4289                        }
     4290                }
     4291               
     4292                /*Add gradje_g_gaussian vector to gradje_g: */
     4293                for( i=0; i<numgrids; i++)grade_g[i]+=grade_g_gaussian[i];
     4294        }
     4295
     4296
     4297        /*Add grade_g to the inputs of this element: */
     4298        this->inputs->AddInput(new TriaVertexInput(GradientEnum,&grade_g[0]));
     4299
     4300        cleanup_and_return:
     4301        xfree((void**)&first_gauss_area_coord);
     4302        xfree((void**)&second_gauss_area_coord);
     4303        xfree((void**)&third_gauss_area_coord);
     4304        xfree((void**)&gauss_weights);
     4305        delete friction;
     4306
     4307}
     4308/*}}}*/
     4309/*FUNCTION Tria::GradjDragStokes {{{1*/
     4310void  Tria::GradjDragStokes(Vec grad_g){
     4311
     4312        int i;
     4313
     4314        /* node data: */
     4315        const int    numgrids=3;
     4316        const int    NDOF2=2;
     4317        double       xyz_list[numgrids][3];
     4318        int          doflist1[numgrids];
     4319        double       dh1dh3[NDOF2][numgrids];
     4320
     4321        /* grid data: */
     4322        double drag;
     4323        double alpha_complement;
     4324        Friction* friction=NULL;
     4325
     4326        /* gaussian points: */
     4327        int     num_gauss,ig;
     4328        double* first_gauss_area_coord  =  NULL;
     4329        double* second_gauss_area_coord =  NULL;
     4330        double* third_gauss_area_coord  =  NULL;
     4331        double* gauss_weights           =  NULL;
     4332        double  gauss_weight;
     4333        double  gauss_l1l2l3[3];
     4334        double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
     4335
     4336        /* parameters: */
     4337        double  vx,vy,vz;
     4338        double  lambda,mu,xi;
     4339        double  bed,thickness,Neff;
     4340        double  surface_normal[3];
     4341        double  bed_normal[3];
     4342        double  dk[NDOF2];
     4343
     4344        /*element vector at the gaussian points: */
     4345        double  grade_g[numgrids]={0.0};
     4346        double  grade_g_gaussian[numgrids];
     4347
     4348        /* Jacobian: */
     4349        double Jdet;
     4350
     4351        /*nodal functions: */
     4352        double l1l2l3[3];
     4353
     4354        /* strain rate: */
     4355        double epsilon[3]; /* epsilon=[exx,eyy,exy];*/
     4356
     4357        /*inputs: */
     4358        bool shelf;
     4359        int  drag_type;
     4360
     4361        /*parameters: */
     4362        double  cm_noisedmp;
     4363        double  cm_mindmp_slope;
     4364        double  cm_mindmp_value;
     4365        double  cm_maxdmp_value;
     4366        double  cm_maxdmp_slope;
     4367
     4368        int analysis_type;
     4369
     4370        /*retrive parameters: */
     4371        parameters->FindParam(&analysis_type,AnalysisTypeEnum);
     4372
     4373        /*retrieve inputs :*/
     4374        inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
     4375        inputs->GetParameterValue(&drag_type,DragTypeEnum);
     4376
     4377        /*retrieve some parameters: */
     4378        this->parameters->FindParam(&cm_noisedmp,CmNoiseDmpEnum);
     4379        this->parameters->FindParam(&cm_mindmp_value,CmMinDmpValueEnum);
     4380        this->parameters->FindParam(&cm_mindmp_slope,CmMinDmpSlopeEnum);
     4381        this->parameters->FindParam(&cm_maxdmp_value,CmMaxDmpValueEnum);
     4382        this->parameters->FindParam(&cm_maxdmp_slope,CmMaxDmpSlopeEnum);
     4383
     4384        /*Get out if shelf*/
     4385        if(shelf)return;
     4386
     4387        /* Get node coordinates and dof list: */
     4388        GetVerticesCoordinates(&xyz_list[0][0], nodes, numgrids);
     4389        GetDofList1(&doflist1[0]);
     4390
     4391        /*Build frictoin element, needed later: */
     4392        inputs->GetParameterValue(&drag_type,DragTypeEnum);
     4393        friction=new Friction("2d",inputs,matpar,analysis_type);
     4394
     4395        /* Get gaussian points and weights (make this a statically initialized list of points? fstd): */
     4396        GaussTria( &num_gauss, &first_gauss_area_coord, &second_gauss_area_coord, &third_gauss_area_coord, &gauss_weights, 4);
     4397
     4398        /* Start  looping on the number of gaussian points: */
     4399        for (ig=0; ig<num_gauss; ig++){
     4400                /*Pick up the gaussian point: */
     4401                gauss_weight=*(gauss_weights+ig);
     4402                gauss_l1l2l3[0]=*(first_gauss_area_coord+ig);
     4403                gauss_l1l2l3[1]=*(second_gauss_area_coord+ig);
     4404                gauss_l1l2l3[2]=*(third_gauss_area_coord+ig);
     4405
     4406                /*Recover alpha_complement and drag: */
     4407                if (drag_type==2) friction->GetAlphaComplement(&alpha_complement, gauss_l1l2l3,VxAverageEnum,VyAverageEnum);
     4408                else alpha_complement=0;
     4409                inputs->GetParameterValue(&drag, &gauss_l1l2l3[0],DragCoefficientEnum);
     4410
     4411                /*recover lambda mu and xi: */
     4412                inputs->GetParameterValue(&lambda, &gauss_l1l2l3[0],AdjointxEnum);
     4413                inputs->GetParameterValue(&mu, &gauss_l1l2l3[0],AdjointyEnum);
     4414                inputs->GetParameterValue(&xi, &gauss_l1l2l3[0],AdjointzEnum);
     4415
     4416                /*recover vx vy and vz: */
     4417                inputs->GetParameterValue(&vx, &gauss_l1l2l3[0],VxEnum);
     4418                inputs->GetParameterValue(&vy, &gauss_l1l2l3[0],VyEnum);
     4419                inputs->GetParameterValue(&vz, &gauss_l1l2l3[0],VzEnum);
     4420
     4421                /*Get normal vector to the bed */
     4422                SurfaceNormal(&surface_normal[0],xyz_list);
     4423
     4424                bed_normal[0]=-surface_normal[0]; //Program is for surface, so the normal to the bed is the opposite of the result
     4425                bed_normal[1]=-surface_normal[1];
     4426                bed_normal[2]=-surface_normal[2];
     4427
     4428                /* Get Jacobian determinant: */
     4429                GetJacobianDeterminant3d(&Jdet, &xyz_list[0][0],gauss_l1l2l3);
     4430
     4431                /* Get nodal functions value at gaussian point:*/
     4432                GetNodalFunctions(l1l2l3, gauss_l1l2l3);
     4433
     4434                /*Get nodal functions derivatives*/
     4435                GetNodalFunctionsDerivatives(&dh1dh3[0][0],&xyz_list[0][0],gauss_l1l2l3);
     4436
     4437                /*Get k derivative: dk/dx */
     4438                inputs->GetParameterDerivativeValue(&dk[0],&xyz_list[0][0],&gauss_l1l2l3[0],DragCoefficientEnum);
     4439
     4440                /*Build gradje_g_gaussian vector (actually -dJ/ddrag): */
     4441                for (i=0;i<numgrids;i++){
     4442                        //standard gradient dJ/dki
     4443                        grade_g_gaussian[i]=(
     4444                                                -lambda*(2*drag*alpha_complement*(vx - vz*bed_normal[0]*bed_normal[2]))
     4445                                                -mu    *(2*drag*alpha_complement*(vy - vz*bed_normal[1]*bed_normal[2]))
     4446                                                -xi    *(2*drag*alpha_complement*(-vx*bed_normal[0]*bed_normal[2]-vy*bed_normal[1]*bed_normal[2]))
     4447                                                )*Jdet*gauss_weight*l1l2l3[i];
     4448
     4449                        //Add regularization term
     4450                        grade_g_gaussian[i]+= - cm_noisedmp*Jdet*gauss_weight*(dh1dh3[0][i]*dk[0]+dh1dh3[1][i]*dk[1]);
     4451
     4452                        //min dampening
     4453                        if(drag<cm_mindmp_value){
     4454                                grade_g_gaussian[i]+= cm_mindmp_slope*Jdet*gauss_weight*l1l2l3[i];
     4455                        }
     4456
     4457                        //max dampening
     4458                        if(drag>cm_maxdmp_value){
     4459                                grade_g_gaussian[i]+= - cm_maxdmp_slope*Jdet*gauss_weight*l1l2l3[i];
     4460                        }
     4461                }
     4462
     4463                /*Add gradje_g_gaussian vector to gradje_g: */
     4464                for( i=0; i<numgrids; i++)grade_g[i]+=grade_g_gaussian[i];
     4465        }
     4466
     4467        /*Add grade_g to the inputs of this element: */
     4468        this->inputs->AddInput(new TriaVertexInput(GradientEnum,&grade_g[0]));
     4469
     4470        cleanup_and_return:
     4471        xfree((void**)&first_gauss_area_coord);
     4472        xfree((void**)&second_gauss_area_coord);
     4473        xfree((void**)&third_gauss_area_coord);
     4474        xfree((void**)&gauss_weights);
     4475        delete friction;
     4476
    41454477}
    41464478/*}}}*/
    … …  
    42814613        }
    42824614
    4283         /*Add grade_g to global vector grad_g: */
    4284         VecSetValues(grad_g,numgrids,doflist1,(const double*)grade_g,ADD_VALUES);
    4285 
    4286 cleanup_and_return:
    4287         xfree((void**)&first_gauss_area_coord);
    4288         xfree((void**)&second_gauss_area_coord);
    4289         xfree((void**)&third_gauss_area_coord);
    4290         xfree((void**)&gauss_weights);
    4291 }
    4292 /*}}}*/
    4293 /*FUNCTION Tria::GradjDrag {{{1*/
    4294 void  Tria::GradjDrag(Vec grad_g){
    4295 
    4296 
    4297         int i;
    4298 
    4299         /* node data: */
    4300         const int    numgrids=3;
    4301         const int    NDOF2=2;
    4302         const int    numdof=NDOF2*numgrids;
    4303         double       xyz_list[numgrids][3];
    4304         int          doflist1[numgrids];
    4305         double       dh1dh3[NDOF2][numgrids];
    4306 
    4307         /* grid data: */
    4308         double adjx_list[numgrids];
    4309         double adjy_list[numgrids];
    4310 
    4311         /* gaussian points: */
    4312         int     num_gauss,ig;
    4313         double* first_gauss_area_coord  =  NULL;
    4314         double* second_gauss_area_coord =  NULL;
    4315         double* third_gauss_area_coord  =  NULL;
    4316         double* gauss_weights           =  NULL;
    4317         double  gauss_weight;
    4318         double  gauss_l1l2l3[3];
    4319 
    4320         /* parameters: */
    4321         double  dk[NDOF2];
    4322         double  vx,vy;
    4323         double  lambda,mu;
    4324         double  bed,thickness,Neff;
    4325         double  alpha_complement;
    4326         int     drag_type;
    4327         double  drag;
    4328         Friction* friction=NULL;
    4329 
    4330         /*element vector at the gaussian points: */
    4331         double  grade_g[numgrids]={0.0};
    4332         double  grade_g_gaussian[numgrids];
    4333 
    4334         /* Jacobian: */
    4335         double Jdet;
    4336 
    4337         /*nodal functions: */
    4338         double l1l2l3[3];
    4339 
    4340         /* strain rate: */
    4341         double epsilon[3]; /* epsilon=[exx,eyy,exy];*/
    4342 
    4343         /*inputs: */
    4344         bool shelf;
    4345 
    4346         /*parameters: */
    4347         double  cm_noisedmp;
    4348         double  cm_mindmp_slope;
    4349         double  cm_mindmp_value;
    4350         double  cm_maxdmp_value;
    4351         double  cm_maxdmp_slope;
    4352 
    4353         int analysis_type;
    4354 
    4355         /*retrive parameters: */
    4356         parameters->FindParam(&analysis_type,AnalysisTypeEnum);
    4357 
    4358         /*retrieve inputs :*/
    4359         inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
    4360 
    4361         /*retrieve some parameters: */
    4362         this->parameters->FindParam(&cm_noisedmp,CmNoiseDmpEnum);
    4363         this->parameters->FindParam(&cm_mindmp_value,CmMinDmpValueEnum);
    4364         this->parameters->FindParam(&cm_mindmp_slope,CmMinDmpSlopeEnum);
    4365         this->parameters->FindParam(&cm_maxdmp_value,CmMaxDmpValueEnum);
    4366         this->parameters->FindParam(&cm_maxdmp_slope,CmMaxDmpSlopeEnum);
    4367 
    4368 
    4369         /*Get out if shelf*/
    4370         if(shelf)return;
    4371 
    4372         /* Get node coordinates and dof list: */
    4373         GetVerticesCoordinates(&xyz_list[0][0], nodes, numgrids);
    4374         GetDofList1(&doflist1[0]);
    4375 
    4376         /*Build frictoin element, needed later: */
    4377         inputs->GetParameterValue(&drag_type,DragTypeEnum);
    4378         friction=new Friction("2d",inputs,matpar,analysis_type);
    4379 
    4380         /* Get gaussian points and weights (make this a statically initialized list of points? fstd): */
    4381         GaussTria( &num_gauss, &first_gauss_area_coord, &second_gauss_area_coord, &third_gauss_area_coord, &gauss_weights, 4);
    4382 
    4383         /* Start  looping on the number of gaussian points: */
    4384         for (ig=0; ig<num_gauss; ig++){
    4385                 /*Pick up the gaussian point: */
    4386                 gauss_weight=*(gauss_weights+ig);
    4387                 gauss_l1l2l3[0]=*(first_gauss_area_coord+ig);
    4388                 gauss_l1l2l3[1]=*(second_gauss_area_coord+ig);
    4389                 gauss_l1l2l3[2]=*(third_gauss_area_coord+ig);
    4390 
    4391                 /*Build alpha_complement_list: */
    4392                 if (drag_type==2) friction->GetAlphaComplement(&alpha_complement, gauss_l1l2l3,VxAverageEnum,VyAverageEnum);
    4393                 else alpha_complement=0;
    4394        
    4395                 /*Recover alpha_complement and k: */
    4396                 inputs->GetParameterValue(&drag, gauss_l1l2l3,DragCoefficientEnum);
    4397 
    4398                 /*recover lambda and mu: */
    4399                 inputs->GetParameterValue(&lambda, gauss_l1l2l3,AdjointxEnum);
    4400                 inputs->GetParameterValue(&mu, gauss_l1l2l3,AdjointyEnum);
    4401                        
    4402                 /*recover vx and vy: */
    4403                 inputs->GetParameterValue(&vx, gauss_l1l2l3,VxEnum);
    4404                 inputs->GetParameterValue(&vy, gauss_l1l2l3,VyEnum);
    4405 
    4406                 /* Get Jacobian determinant: */
    4407                 GetJacobianDeterminant2d(&Jdet, &xyz_list[0][0],gauss_l1l2l3);
    4408                
    4409                 /* Get nodal functions value at gaussian point:*/
    4410                 GetNodalFunctions(l1l2l3, gauss_l1l2l3);
    4411 
    4412                 /*Get nodal functions derivatives*/
    4413                 GetNodalFunctionsDerivatives(&dh1dh3[0][0],&xyz_list[0][0],gauss_l1l2l3);
    4414 
    4415                 /*Get k derivative: dk/dx */
    4416                 inputs->GetParameterDerivativeValue(&dk[0],&xyz_list[0][0],&gauss_l1l2l3[0],DragCoefficientEnum);
    4417 
    4418                 /*Build gradje_g_gaussian vector (actually -dJ/ddrag): */
    4419                 for (i=0;i<numgrids;i++){
    4420 
    4421                         //standard term dJ/dki
    4422                         grade_g_gaussian[i]=-2*drag*alpha_complement*((lambda*vx+mu*vy))*Jdet*gauss_weight*l1l2l3[i];
    4423 
    4424                         //noise dampening d/dki(1/2*(dk/dx)^2)
    4425                         grade_g_gaussian[i]+=-cm_noisedmp*Jdet*gauss_weight*(dh1dh3[0][i]*dk[0]+dh1dh3[1][i]*dk[1]);
    4426                        
    4427                         //min dampening
    4428                         if(drag<cm_mindmp_value){
    4429                                 grade_g_gaussian[i]+=cm_mindmp_slope*Jdet*gauss_weight*l1l2l3[i];
    4430                         }
    4431 
    4432                         //max dampening
    4433                         if(drag>cm_maxdmp_value){
    4434                                 grade_g_gaussian[i]+= - cm_maxdmp_slope*Jdet*gauss_weight*l1l2l3[i];
    4435                         }
    4436                 }
    4437                
    4438                 /*Add gradje_g_gaussian vector to gradje_g: */
    4439                 for( i=0; i<numgrids; i++)grade_g[i]+=grade_g_gaussian[i];
    4440         }
    4441 
    4442         /*Add grade_g to global vector grad_g: */
    4443         VecSetValues(grad_g,numgrids,doflist1,(const double*)grade_g,ADD_VALUES);
    4444        
     4615        /*Add grade_g to the inputs of this element: */
     4616        this->inputs->AddInput(new TriaVertexInput(GradientEnum,&grade_g[0]));
     4617
    44454618        cleanup_and_return:
    44464619        xfree((void**)&first_gauss_area_coord);
    … …  
    44484621        xfree((void**)&third_gauss_area_coord);
    44494622        xfree((void**)&gauss_weights);
    4450         delete friction;
    4451 
    4452 }
    4453 /*}}}*/
    4454 /*FUNCTION Tria::GradjDragStokes {{{1*/
    4455 void  Tria::GradjDragStokes(Vec grad_g){
    4456 
    4457         int i;
    4458 
    4459         /* node data: */
    4460         const int    numgrids=3;
    4461         const int    NDOF2=2;
    4462         double       xyz_list[numgrids][3];
    4463         int          doflist1[numgrids];
    4464         double       dh1dh3[NDOF2][numgrids];
    4465 
    4466         /* grid data: */
    4467         double drag;
    4468         double alpha_complement;
    4469         Friction* friction=NULL;
    4470 
    4471         /* gaussian points: */
    4472         int     num_gauss,ig;
    4473         double* first_gauss_area_coord  =  NULL;
    4474         double* second_gauss_area_coord =  NULL;
    4475         double* third_gauss_area_coord  =  NULL;
    4476         double* gauss_weights           =  NULL;
    4477         double  gauss_weight;
    4478         double  gauss_l1l2l3[3];
    4479         double  gaussgrids[numgrids][numgrids]={{1,0,0},{0,1,0},{0,0,1}};
    4480 
    4481         /* parameters: */
    4482         double  vx,vy,vz;
    4483         double  lambda,mu,xi;
    4484         double  bed,thickness,Neff;
    4485         double  surface_normal[3];
    4486         double  bed_normal[3];
    4487         double  dk[NDOF2];
    4488 
    4489         /*element vector at the gaussian points: */
    4490         double  grade_g[numgrids]={0.0};
    4491         double  grade_g_gaussian[numgrids];
    4492 
    4493         /* Jacobian: */
    4494         double Jdet;
    4495 
    4496         /*nodal functions: */
    4497         double l1l2l3[3];
    4498 
    4499         /* strain rate: */
    4500         double epsilon[3]; /* epsilon=[exx,eyy,exy];*/
    4501 
    4502         /*inputs: */
    4503         bool shelf;
    4504         int  drag_type;
    4505 
    4506         /*parameters: */
    4507         double  cm_noisedmp;
    4508         double  cm_mindmp_slope;
    4509         double  cm_mindmp_value;
    4510         double  cm_maxdmp_value;
    4511         double  cm_maxdmp_slope;
    4512 
    4513         int analysis_type;
    4514 
    4515         /*retrive parameters: */
    4516         parameters->FindParam(&analysis_type,AnalysisTypeEnum);
    4517 
    4518         /*retrieve inputs :*/
    4519         inputs->GetParameterValue(&shelf,ElementOnIceShelfEnum);
    4520         inputs->GetParameterValue(&drag_type,DragTypeEnum);
    4521 
    4522         /*retrieve some parameters: */
    4523         this->parameters->FindParam(&cm_noisedmp,CmNoiseDmpEnum);
    4524         this->parameters->FindParam(&cm_mindmp_value,CmMinDmpValueEnum);
    4525         this->parameters->FindParam(&cm_mindmp_slope,CmMinDmpSlopeEnum);
    4526         this->parameters->FindParam(&cm_maxdmp_value,CmMaxDmpValueEnum);
    4527         this->parameters->FindParam(&cm_maxdmp_slope,CmMaxDmpSlopeEnum);
    4528 
    4529         /*Get out if shelf*/
    4530         if(shelf)return;
    4531 
    4532         /* Get node coordinates and dof list: */
    4533         GetVerticesCoordinates(&xyz_list[0][0], nodes, numgrids);
    4534         GetDofList1(&doflist1[0]);
    4535 
    4536         /*Build frictoin element, needed later: */
    4537         inputs->GetParameterValue(&drag_type,DragTypeEnum);
    4538         friction=new Friction("2d",inputs,matpar,analysis_type);
    4539 
    4540         /* Get gaussian points and weights (make this a statically initialized list of points? fstd): */
    4541         GaussTria( &num_gauss, &first_gauss_area_coord, &second_gauss_area_coord, &third_gauss_area_coord, &gauss_weights, 4);
    4542 
    4543         /* Start  looping on the number of gaussian points: */
    4544         for (ig=0; ig<num_gauss; ig++){
    4545                 /*Pick up the gaussian point: */
    4546                 gauss_weight=*(gauss_weights+ig);
    4547                 gauss_l1l2l3[0]=*(first_gauss_area_coord+ig);
    4548                 gauss_l1l2l3[1]=*(second_gauss_area_coord+ig);
    4549                 gauss_l1l2l3[2]=*(third_gauss_area_coord+ig);
    4550 
    4551                 /*Recover alpha_complement and drag: */
    4552                 if (drag_type==2) friction->GetAlphaComplement(&alpha_complement, gauss_l1l2l3,VxAverageEnum,VyAverageEnum);
    4553                 else alpha_complement=0;
    4554                 inputs->GetParameterValue(&drag, &gauss_l1l2l3[0],DragCoefficientEnum);
    4555 
    4556                 /*recover lambda mu and xi: */
    4557                 inputs->GetParameterValue(&lambda, &gauss_l1l2l3[0],AdjointxEnum);
    4558                 inputs->GetParameterValue(&mu, &gauss_l1l2l3[0],AdjointyEnum);
    4559                 inputs->GetParameterValue(&xi, &gauss_l1l2l3[0],AdjointzEnum);
    4560 
    4561                 /*recover vx vy and vz: */
    4562                 inputs->GetParameterValue(&vx, &gauss_l1l2l3[0],VxEnum);
    4563                 inputs->GetParameterValue(&vy, &gauss_l1l2l3[0],VyEnum);
    4564                 inputs->GetParameterValue(&vz, &gauss_l1l2l3[0],VzEnum);
    4565 
    4566                 /*Get normal vector to the bed */
    4567                 SurfaceNormal(&surface_normal[0],xyz_list);
    4568 
    4569                 bed_normal[0]=-surface_normal[0]; //Program is for surface, so the normal to the bed is the opposite of the result
    4570                 bed_normal[1]=-surface_normal[1];
    4571                 bed_normal[2]=-surface_normal[2];
    4572 
    4573                 /* Get Jacobian determinant: */
    4574                 GetJacobianDeterminant3d(&Jdet, &xyz_list[0][0],gauss_l1l2l3);
    4575 
    4576                 /* Get nodal functions value at gaussian point:*/
    4577                 GetNodalFunctions(l1l2l3, gauss_l1l2l3);
    4578 
    4579                 /*Get nodal functions derivatives*/
    4580                 GetNodalFunctionsDerivatives(&dh1dh3[0][0],&xyz_list[0][0],gauss_l1l2l3);
    4581 
    4582                 /*Get k derivative: dk/dx */
    4583                 inputs->GetParameterDerivativeValue(&dk[0],&xyz_list[0][0],&gauss_l1l2l3[0],DragCoefficientEnum);
    4584 
    4585                 /*Build gradje_g_gaussian vector (actually -dJ/ddrag): */
    4586                 for (i=0;i<numgrids;i++){
    4587                         //standard gradient dJ/dki
    4588                         grade_g_gaussian[i]=(
    4589                                                 -lambda*(2*drag*alpha_complement*(vx - vz*bed_normal[0]*bed_normal[2]))
    4590                                                 -mu    *(2*drag*alpha_complement*(vy - vz*bed_normal[1]*bed_normal[2]))
    4591                                                 -xi    *(2*drag*alpha_complement*(-vx*bed_normal[0]*bed_normal[2]-vy*bed_normal[1]*bed_normal[2]))
    4592                                                 )*Jdet*gauss_weight*l1l2l3[i];
    4593 
    4594                         //Add regularization term
    4595                         grade_g_gaussian[i]+= - cm_noisedmp*Jdet*gauss_weight*(dh1dh3[0][i]*dk[0]+dh1dh3[1][i]*dk[1]);
    4596 
    4597                         //min dampening
    4598                         if(drag<cm_mindmp_value){
    4599                                 grade_g_gaussian[i]+= cm_mindmp_slope*Jdet*gauss_weight*l1l2l3[i];
    4600                         }
    4601 
    4602                         //max dampening
    4603                         if(drag>cm_maxdmp_value){
    4604                                 grade_g_gaussian[i]+= - cm_maxdmp_slope*Jdet*gauss_weight*l1l2l3[i];
    4605                         }
    4606                 }
    4607 
    4608                 /*Add gradje_g_gaussian vector to gradje_g: */
    4609                 for( i=0; i<numgrids; i++)grade_g[i]+=grade_g_gaussian[i];
    4610         }
    4611 
    4612         /*Add grade_g to global vector grad_g: */
    4613         VecSetValues(grad_g,numgrids,doflist1,(const double*)grade_g,ADD_VALUES);
    4614 
    4615         cleanup_and_return:
    4616         xfree((void**)&first_gauss_area_coord);
    4617         xfree((void**)&second_gauss_area_coord);
    4618         xfree((void**)&third_gauss_area_coord);
    4619         xfree((void**)&gauss_weights);
    4620         delete friction;
    4621 
    46224623}
    46234624/*}}}*/
    … …  
    54375438}
    54385439/*}}}*/
     5440/*FUNCTION Tria::ScaleInput(int enum_type,double scale_factor){{{1*/
     5441void  Tria::ScaleInput(int enum_type,double scale_factor){
     5442
     5443        Input* input=NULL;
     5444
     5445        /*Make a copy of the original input: */
     5446        input=(Input*)this->inputs->GetInput(enum_type);
     5447
     5448        /*Scale: */
     5449        input->Scale(scale_factor);
     5450}
     5451/*}}}*/

  • issm/trunk/src/c/objects/Elements/Tria.h

    r4046 r4047  
    8888                void  GetJacobian(double* J, double* xyz_list,double* gauss_l1l2l3);
    8989                void  Du(Vec du_g);
    90                 void  Gradj(Vec grad_g,int control_type);
    91                 void  GradjDrag(Vec grad_g);
    92                 void  GradjDragStokes(Vec grad_g);
     90                void  Gradj(int control_type);
     91                void  GradjDrag(void);
     92                void  GradjDragStokes(void);
     93                void  GradjB(void);
    9394                void  SurfaceNormal(double* surface_normal, double xyz_list[3][3]);
    94                 void  GradjB(Vec grad_g);
    9595                double Misfit(void);
    9696                double SurfaceArea(void);
    … …  
    140140                void  MaxAbsVz(double* pmaxabsvz, bool process_units);
    141141                void  DuplicateInput(int original_enum,int new_enum);
     142                void  ScaleInput(int enum_type,double scale_factor);
    142143
    143144

  • issm/trunk/src/c/objects/Inputs/BeamVertexInput.cpp

    r4043 r4047  
    250250}
    251251/*}}}*/
     252/*FUNCTION BeamVertexInput::Scale(double scale_factor){{{1*/
     253void BeamVertexInput::Scale(double scale_factor){
     254
     255        int i;
     256        const int numgrids=2;
     257
     258        for(i=0;i<numgrids;i++)values[i]=values[i]*scale_factor;
     259}
     260/*}}}*/

  • issm/trunk/src/c/objects/Inputs/BeamVertexInput.h

    r4043 r4047  
    7474                void ChangeEnum(int newenumtype);
    7575                void SquareMin(double* psquaremin, bool process_units,Parameters* parameters);
     76                void Scale(double scale_factor);
    7677                /*}}}*/
    7778

  • issm/trunk/src/c/objects/Inputs/BoolInput.cpp

    r4042 r4047  
    223223}
    224224/*}}}*/
     225/*FUNCTION BoolInput::Scale(double scale_factor){{{1*/
     226void BoolInput::Scale(double scale_factor){
     227        /*a bool cannot be scaled: */
     228}
     229/*}}}*/

  • issm/trunk/src/c/objects/Inputs/BoolInput.h

    r4043 r4047  
    7474                void ChangeEnum(int newenumtype);
    7575                void SquareMin(double* psquaremin, bool process_units,Parameters* parameters);
     76                void Scale(double scale_factor);
    7677                /*}}}*/
    7778

  • issm/trunk/src/c/objects/Inputs/DoubleInput.cpp

    r4042 r4047  
    234234}
    235235/*}}}*/
     236/*FUNCTION DoubleInput::Scale(double scale_factor){{{1*/
     237void DoubleInput::Scale(double scale_factor){
     238        value=value*scale_factor;
     239}
     240/*}}}*/

  • issm/trunk/src/c/objects/Inputs/DoubleInput.h

    r4043 r4047  
    7575                void ChangeEnum(int newenumtype);
    7676                void SquareMin(double* psquaremin, bool process_units,Parameters* parameters);
     77                void Scale(double scale_factor);
    7778                /*}}}*/
    7879

  • issm/trunk/src/c/objects/Inputs/Input.h

    r4042 r4047  
    5151                virtual Result* SpawnResult(int step, double time)=0;
    5252                virtual void   SquareMin(double* psquaremin, bool process_units,Parameters* parameters)=0;
     53                virtual void   Scale(double scale_factor)=0;
    5354                /*}}}*/
    5455

  • issm/trunk/src/c/objects/Inputs/IntInput.cpp

    r4042 r4047  
    221221}
    222222/*}}}*/
     223/*FUNCTION IntInput::Scale(double scale_factor){{{1*/
     224void IntInput::Scale(double scale_factor){
     225        value=value*scale_factor;
     226}
     227/*}}}*/

  • issm/trunk/src/c/objects/Inputs/IntInput.h

    r4043 r4047  
    7474                void ChangeEnum(int newenumtype);
    7575                void SquareMin(double* psquaremin, bool process_units,Parameters* parameters);
     76                void Scale(double scale_factor);
    7677                /*}}}*/
    7778

  • issm/trunk/src/c/objects/Inputs/PentaVertexInput.cpp

    r4043 r4047  
    899899}
    900900/*}}}*/
     901/*FUNCTION PentaVertexInput::Scale(double scale_factor){{{1*/
     902void PentaVertexInput::Scale(double scale_factor){
     903       
     904        int i;
     905        const int numgrids=6;
     906
     907        for(i=0;i<numgrids;i++)values[i]=values[i]*scale_factor;
     908}
     909/*}}}*/

  • issm/trunk/src/c/objects/Inputs/PentaVertexInput.h

    r4043 r4047  
    8383                void GetBStokes(double* B, double* xyz_list, double* gauss_coord);
    8484                void SquareMin(double* psquaremin, bool process_units,Parameters* parameters);
     85                void Scale(double scale_factor);
    8586                /*}}}*/
    8687

  • issm/trunk/src/c/objects/Inputs/SingVertexInput.cpp

    r4043 r4047  
    224224}
    225225/*}}}*/
     226/*FUNCTION SingVertexInput::Scale(double scale_factor){{{1*/
     227void SingVertexInput::Scale(double scale_factor){
     228        value=value*scale_factor;
     229}
     230/*}}}*/

  • issm/trunk/src/c/objects/Inputs/SingVertexInput.h

    r4043 r4047  
    7373                void ChangeEnum(int newenumtype);
    7474                void SquareMin(double* psquaremin, bool process_units,Parameters* parameters);
     75                void Scale(double scale_factor);
    7576                /*}}}*/
    7677

  • issm/trunk/src/c/objects/Inputs/TriaVertexInput.cpp

    r4043 r4047  
    473473}
    474474/*}}}*/
     475/*FUNCTION TriaVertexInput::Scale(double scale_factor){{{1*/
     476void TriaVertexInput::Scale(double scale_factor){
     477       
     478        int i;
     479        const int numgrids=3;
     480
     481        for(i=0;i<numgrids;i++)values[i]=values[i]*scale_factor;
     482}
     483/*}}}*/

  • issm/trunk/src/c/objects/Inputs/TriaVertexInput.h

    r4043 r4047  
    8080                void GetJacobianInvert(double*  Jinv, double* xyz_list,double* gauss);
    8181                void SquareMin(double* psquaremin, bool process_units,Parameters* parameters);
     82                void Scale(double scale_factor);
    8283                /*}}}*/
    8384

  • issm/trunk/src/c/solutions/control_core.cpp

    r4043 r4047  
    33 */
    44
     5#include "../toolkits/toolkits.h"
     6#include "../objects/objects.h"
     7#include "../shared/shared.h"
     8#include "../EnumDefinitions/EnumDefinitions.h"
    59#include "./solutions.h"
    610#include "../modules/modules.h"
    7 #include "../EnumDefinitions/EnumDefinitions.h"
     11#include "../include/include.h"
    812#include "../solvers/solvers.h"
    913
    10 Results* control_core(Model* model){
    11 
    12         extern int my_rank;
    13 
    14         /*fem model: */
    15         FemModel* fem_model=NULL;
    16 
    17         /*output: */
    18         Results* results=NULL;
    19 
    20         /*Intermediary: */
    21         Results* diagnostic_results=NULL;
    22         Results* gradjcompute_results=NULL;
    23         Results* steadystate_results=NULL;
     14void control_core(FemModel* femmodel){
     15
    2416        Vec     u_g=NULL;
    2517        Vec     t_g=NULL;
    … …  
    2921        double* fit=NULL;
    3022        double* optscal=NULL;
    31         int     gsize;
    3223        double* maxiter=NULL;
    3324        double* cm_jump=NULL;
    … …  
    5647        int numberofnodes;
    5748
    58         //initialize results
    59         results=new Results();
    60 
    6149        /*Process models*/
    62         ControlInitialization(model);
    63         fem_model=model->GetActiveFormulation();
     50        ControlInitialization(femmodel);
    6451
    6552        /*Recover parameters used throughout the solution:*/
    66         model->FindParam(&nsteps,NStepsEnum);
    67         model->FindParam(&control_type,ControlTypeEnum);
    68         model->FindParam(&fit,NULL,NULL,FitEnum);
    69         model->FindParam(&optscal,NULL,NULL,OptScalEnum);
    70         model->FindParam(&maxiter,NULL,NULL,MaxIterEnum);
    71         model->FindParam(&cm_jump,NULL,NULL,CmJumpEnum);
    72         model->FindParam(&eps_cm,EpsCmEnum);
    73         model->FindParam(&tolx,TolXEnum);
    74         model->FindParam(&cm_min,CmMinEnum);
    75         model->FindParam(&cm_max,CmMaxEnum);
    76         model->FindParam(&cm_gradient,CmGradientEnum);
    77         model->FindParam(&param_g,NULL,NULL,ControlParameterEnum);
    78         model->FindParam(&analysis_type,AnalysisTypeEnum);
    79         model->FindParam(&sub_analysis_type,SubAnalysisTypeEnum);
    80         model->FindParam(&numberofnodes,NumberOfNodesEnum);
    81         model->FindParam(&control_steady,ControlSteadyEnum);
    82         gsize=fem_model->nodes->NumberOfDofs();
     53        femmodel->parameters->FindParam(&nsteps,NStepsEnum);
     54        femmodel->parameters->FindParam(&control_type,ControlTypeEnum);
     55        femmodel->parameters->FindParam(&fit,NULL,NULL,FitEnum);
     56        femmodel->parameters->FindParam(&optscal,NULL,NULL,OptScalEnum);
     57        femmodel->parameters->FindParam(&maxiter,NULL,NULL,MaxIterEnum);
     58        femmodel->parameters->FindParam(&cm_jump,NULL,NULL,CmJumpEnum);
     59        femmodel->parameters->FindParam(&eps_cm,EpsCmEnum);
     60        femmodel->parameters->FindParam(&tolx,TolXEnum);
     61        femmodel->parameters->FindParam(&cm_min,CmMinEnum);
     62        femmodel->parameters->FindParam(&cm_max,CmMaxEnum);
     63        femmodel->parameters->FindParam(&cm_gradient,CmGradientEnum);
     64        femmodel->parameters->FindParam(&param_g,NULL,NULL,ControlParameterEnum);
     65        femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum);
     66        femmodel->parameters->FindParam(&sub_analysis_type,SubAnalysisTypeEnum);
     67        femmodel->parameters->FindParam(&numberofnodes,NumberOfNodesEnum);
     68        femmodel->parameters->FindParam(&control_steady,ControlSteadyEnum);
    8369
    8470        /*Initialize misfit: */
    … …  
    8975
    9076                _printf_("\n%s%i%s%i\n","   control method step ",n+1,"/",nsteps);
    91                 model->UpdateInputsFromVector(param_g,control_type,VertexEnum);
    92                 model->UpdateInputsFromConstant(fit[n],FitEnum);
     77                inputs->Add(control_type,param_g,1,numberofnodes);
     78                femmodel->UpdateInputsFromConstant(fit[n],FitEnum);
    9379               
    9480                /*In case we are running a steady state control method, compute new temperature field using new parameter
    9581                 * distribution: */
    96                 if (control_steady){
    97                         steadystate_results= steadystate_core(model);
    98                         VecFree(&t_g); steadystate_results->FindResult(&t_g,"t_g");
    99                         delete steadystate_results;
    100                         model->UpdateInputsFromVector(t_g,TemperatureEnum,VertexEnum);
    101                 }
     82                if (control_steady) steadystate_core(model);
    10283       
    103                 _printf_("%s\n","      computing gradJ...");
    104                 gradjcompute_results= gradjcompute_core(model);
    105                 gradjcompute_results->FindResult(&grad_g,"grad_g");
    106                 delete gradjcompute_results;
     84                if(verbose)_printf_("%s\n","      computing gradJ...");
     85                gradient_core(femmodel);
    10786
    10887                /*Return gradient if asked: */
    10988                if (cm_gradient){
    110                        
    111                         /*Plug results into output dataset: */
    112                         results->AddObject(new Result(results->Size()+1,0,1,"grad_g",grad_g));
    113                         results->AddObject(new Result(results->Size()+1,0,1,"analysis_type",EnumAsString(analysis_type)));
    114        
    115                         /*Free ressources: */
    116                         xfree((void**)&control_type);
    117                         xfree((void**)&fit);
    118                         xfree((void**)&optscal);
    119                         xfree((void**)&maxiter);
    120                         xfree((void**)&cm_jump);
    121                         xfree((void**)&grad_g_double);
    122                         xfree((void**)&param_g);
    123                         VecFree(&u_g);
    124                         VecFree(&t_g);
    125                         VecFree(&m_g);
    126                         xfree((void**)&J);
    127                         return results;
     89                        /*Transfer gradient from input to results: */
     90                        InputToResultx(femodel->elements,femodel->nodes,femodel->vertices,femodel->loads,femodel->materials,femodel->parameters,GradientEnum);
     91                        return;
    12892                }
    12993

  • issm/trunk/src/c/solutions/diagnostic_core.cpp

    r4046 r4047  
    6868                        //"recondition" pressure computed previously:
    6969                        DuplicateInputx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,PressureEnum,PressureStokesEnum);
    70                         ScaleInputx(femmmodel,PressureStokesEnum,1.0/stokesreconditioning);
     70                        ScaleInputx(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,PressureStokesEnum,1.0/stokesreconditioning);
    7171
    7272                        if(verbose)_printf_("%s\n"," update boundary conditions for stokes using velocities previously computed...");
    … …  
    7474
    7575                        if(verbose)_printf_("%s\n"," computing stokes velocities and pressure ...");
    76                         solver_diagnostic_nonlinear(NULL,NULL,NULL,NULL,fem_ds,DiagnosticAnalysisEnum,StokesAnalysisEnum);
     76                        solver_diagnostic_nonlinear(NULL,NULL,NULL,NULL,femmodel,DiagnosticAnalysisEnum,StokesAnalysisEnum);
    7777                }
    7878        }

  • issm/trunk/src/c/solutions/solutions.h

    r4043 r4047  
    1313
    1414/*cores: */
    15 Results* gradjcompute_core(FemModel* model);
    16 void diagnostic_core(FemModel* model);
    1715Results* prognostic_core(FemModel* model);
    1816Results* prognostic2_core(FemModel* model);
    … …  
    2624Results* transient_core_2d(FemModel* model);
    2725Results* transient_core_3d(FemModel* model);
    28 Results* thermal_core(FemModel* model);
     26void adjoint_core(FemModel* model);
     27void gradient_core(FemModel* model);
     28void diagnostic_core(FemModel* model);
     29void thermal_core(FemModel* model);
    2930void surfaceslope_core(FemModel* femmodel);
    3031void bedslope_core(FemModel* femmodel);

  • issm/trunk/src/c/solvers/solvers.h

    r4029 r4047  
    1515void solver_diagnostic_nonlinear(Vec* pug,Mat* pK_ff0,Mat* pK_fs0, FemModel* fem,bool conserve_loads,int analysis_type,int sub_analysis_type);
    1616void solver_linear(Vec* pug, FemModel* femmodel,int  analysis_type,int sub_analysis_type);
     17void solver_adjoint(Vec* pug, FemModel* femmodel,int  analysis_type,int sub_analysis_type);
    1718
    1819

  • issm/trunk/src/m/enum/ConnectivityEnum.m

    r4043 r4047  
    99%      macro=ConnectivityEnum()
    1010
    11 macro=201;
     11macro=202;

  • issm/trunk/src/m/enum/ControlParameterEnum.m

    r4043 r4047  
    99%      macro=ControlParameterEnum()
    1010
    11 macro=202;
     11macro=203;

  • issm/trunk/src/m/enum/ControlSteadyEnum.m

    r4043 r4047  
    99%      macro=ControlSteadyEnum()
    1010
    11 macro=203;
     11macro=204;

  • issm/trunk/src/m/enum/DakotaParameterEnum.m

    r4043 r4047  
    99%      macro=DakotaParameterEnum()
    1010
    11 macro=204;
     11macro=205;

  • issm/trunk/src/m/enum/DimEnum.m

    r4043 r4047  
    99%      macro=DimEnum()
    1010
    11 macro=205;
     11macro=206;

  • issm/trunk/src/m/enum/EnumAsString.m

    r4043 r4047  
    168168        case PressureOldEnum(), string='PressureOld'; return
    169169        case QmuPressureEnum(), string='QmuPressure'; return
    170         case StokesPressureEnum(), string='StokesPressure'; return
     170        case PressureStokesEnum(), string='PressureStokes'; return
    171171        case ResetPenaltiesEnum(), string='ResetPenalties'; return
    172172        case RheologyBEnum(), string='RheologyB'; return
    … …  
    212212        case CmMaxEnum(), string='CmMax'; return
    213213        case CmMinEnum(), string='CmMin'; return
     214        case GradientEnum(), string='Gradient'; return
    214215        case ConnectivityEnum(), string='Connectivity'; return
    215216        case ControlParameterEnum(), string='ControlParameter'; return

  • issm/trunk/src/m/enum/EpsAbsEnum.m

    r4043 r4047  
    99%      macro=EpsAbsEnum()
    1010
    11 macro=206;
     11macro=207;

  • issm/trunk/src/m/enum/EpsCmEnum.m

    r4043 r4047  
    99%      macro=EpsCmEnum()
    1010
    11 macro=207;
     11macro=208;

  • issm/trunk/src/m/enum/EpsRelEnum.m

    r4043 r4047  
    99%      macro=EpsRelEnum()
    1010
    11 macro=208;
     11macro=209;

  • issm/trunk/src/m/enum/EpsResEnum.m

    r4043 r4047  
    99%      macro=EpsResEnum()
    1010
    11 macro=209;
     11macro=210;

  • issm/trunk/src/m/enum/ExtrudeParamEnum.m

    r4043 r4047  
    99%      macro=ExtrudeParamEnum()
    1010
    11 macro=210;
     11macro=211;

  • issm/trunk/src/m/enum/HeatCapacityEnum.m

    r4043 r4047  
    99%      macro=HeatCapacityEnum()
    1010
    11 macro=211;
     11macro=212;

  • issm/trunk/src/m/enum/IsHutterEnum.m

    r4043 r4047  
    99%      macro=IsHutterEnum()
    1010
    11 macro=212;
     11macro=213;

  • issm/trunk/src/m/enum/IsMacAyealPattynEnum.m

    r4043 r4047  
    99%      macro=IsMacAyealPattynEnum()
    1010
    11 macro=213;
     11macro=214;

  • issm/trunk/src/m/enum/IsStokesEnum.m

    r4043 r4047  
    99%      macro=IsStokesEnum()
    1010
    11 macro=214;
     11macro=215;

  • issm/trunk/src/m/enum/LatentHeatEnum.m

    r4043 r4047  
    99%      macro=LatentHeatEnum()
    1010
    11 macro=215;
     11macro=216;

  • issm/trunk/src/m/enum/LowmemEnum.m

    r4043 r4047  
    99%      macro=LowmemEnum()
    1010
    11 macro=216;
     11macro=217;

  • issm/trunk/src/m/enum/MaxIterEnum.m

    r4043 r4047  
    99%      macro=MaxIterEnum()
    1010
    11 macro=217;
     11macro=218;

  • issm/trunk/src/m/enum/MaxNonlinearIterationsEnum.m

    r4043 r4047  
    99%      macro=MaxNonlinearIterationsEnum()
    1010
    11 macro=218;
     11macro=219;

  • issm/trunk/src/m/enum/MeltingPointEnum.m

    r4043 r4047  
    99%      macro=MeltingPointEnum()
    1010
    11 macro=219;
     11macro=220;

  • issm/trunk/src/m/enum/MinMechanicalConstraintsEnum.m

    r4043 r4047  
    99%      macro=MinMechanicalConstraintsEnum()
    1010
    11 macro=220;
     11macro=221;

  • issm/trunk/src/m/enum/MinThermalConstraintsEnum.m

    r4043 r4047  
    99%      macro=MinThermalConstraintsEnum()
    1010
    11 macro=221;
     11macro=222;

  • issm/trunk/src/m/enum/NStepsEnum.m

    r4043 r4047  
    99%      macro=NStepsEnum()
    1010
    11 macro=222;
     11macro=223;

  • issm/trunk/src/m/enum/NdtEnum.m

    r4043 r4047  
    99%      macro=NdtEnum()
    1010
    11 macro=223;
     11macro=224;

  • issm/trunk/src/m/enum/NumOutputEnum.m

    r4043 r4047  
    99%      macro=NumOutputEnum()
    1010
    11 macro=224;
     11macro=225;

  • issm/trunk/src/m/enum/NumRiftsEnum.m

    r4043 r4047  
    99%      macro=NumRiftsEnum()
    1010
    11 macro=225;
     11macro=226;

  • issm/trunk/src/m/enum/NumberOfDofsPerNodeEnum.m

    r4043 r4047  
    99%      macro=NumberOfDofsPerNodeEnum()
    1010
    11 macro=226;
     11macro=227;

  • issm/trunk/src/m/enum/NumberOfElementsEnum.m

    r4043 r4047  
    99%      macro=NumberOfElementsEnum()
    1010
    11 macro=227;
     11macro=228;

  • issm/trunk/src/m/enum/NumberOfNodesEnum.m

    r4043 r4047  
    99%      macro=NumberOfNodesEnum()
    1010
    11 macro=228;
     11macro=229;

  • issm/trunk/src/m/enum/NumberOfVerticesEnum.m

    r4043 r4047  
    99%      macro=NumberOfVerticesEnum()
    1010
    11 macro=229;
     11macro=230;

  • issm/trunk/src/m/enum/OptScalEnum.m

    r4043 r4047  
    99%      macro=OptScalEnum()
    1010
    11 macro=230;
     11macro=231;

  • issm/trunk/src/m/enum/OutputFileNameEnum.m

    r4043 r4047  
    99%      macro=OutputFileNameEnum()
    1010
    11 macro=231;
     11macro=232;

  • issm/trunk/src/m/enum/ParameterOutputEnum.m

    r4043 r4047  
    99%      macro=ParameterOutputEnum()
    1010
    11 macro=232;
     11macro=233;

  • issm/trunk/src/m/enum/PenaltyMeltingEnum.m

    r4043 r4047  
    99%      macro=PenaltyMeltingEnum()
    1010
    11 macro=233;
     11macro=234;

  • issm/trunk/src/m/enum/QmuAnalysisEnum.m

    r4043 r4047  
    99%      macro=QmuAnalysisEnum()
    1010
    11 macro=234;
     11macro=235;

  • issm/trunk/src/m/enum/QmuErrNameEnum.m

    r4043 r4047  
    99%      macro=QmuErrNameEnum()
    1010
    11 macro=235;
     11macro=236;

  • issm/trunk/src/m/enum/QmuInNameEnum.m

    r4043 r4047  
    99%      macro=QmuInNameEnum()
    1010
    11 macro=236;
     11macro=237;

  • issm/trunk/src/m/enum/QmuMassFluxSegmentsEnum.m

    r4043 r4047  
    99%      macro=QmuMassFluxSegmentsEnum()
    1010
    11 macro=237;
     11macro=238;

  • issm/trunk/src/m/enum/QmuNPartEnum.m

    r4043 r4047  
    99%      macro=QmuNPartEnum()
    1010
    11 macro=238;
     11macro=239;

  • issm/trunk/src/m/enum/QmuOutNameEnum.m

    r4043 r4047  
    99%      macro=QmuOutNameEnum()
    1010
    11 macro=239;
     11macro=240;

  • issm/trunk/src/m/enum/QmuPartEnum.m

    r4043 r4047  
    99%      macro=QmuPartEnum()
    1010
    11 macro=240;
     11macro=241;

  • issm/trunk/src/m/enum/ResponseDescriptorsEnum.m

    r4043 r4047  
    99%      macro=ResponseDescriptorsEnum()
    1010
    11 macro=241;
     11macro=242;

  • issm/trunk/src/m/enum/SolverStringEnum.m

    r4043 r4047  
    99%      macro=SolverStringEnum()
    1010
    11 macro=242;
     11macro=243;

  • issm/trunk/src/m/enum/SparsityEnum.m

    r4043 r4047  
    99%      macro=SparsityEnum()
    1010
    11 macro=243;
     11macro=244;

  • issm/trunk/src/m/enum/StringAsEnum.m

    r4043 r4047  
    166166elseif (strcmpi(name,'PressureOld')), enum=PressureOldEnum(); return
    167167elseif (strcmpi(name,'QmuPressure')), enum=QmuPressureEnum(); return
    168 elseif (strcmpi(name,'StokesPressure')), enum=StokesPressureEnum(); return
     168elseif (strcmpi(name,'PressureStokes')), enum=PressureStokesEnum(); return
    169169elseif (strcmpi(name,'ResetPenalties')), enum=ResetPenaltiesEnum(); return
    170170elseif (strcmpi(name,'RheologyB')), enum=RheologyBEnum(); return
    … …  
    210210elseif (strcmpi(name,'CmMax')), enum=CmMaxEnum(); return
    211211elseif (strcmpi(name,'CmMin')), enum=CmMinEnum(); return
     212elseif (strcmpi(name,'Gradient')), enum=GradientEnum(); return
    212213elseif (strcmpi(name,'Connectivity')), enum=ConnectivityEnum(); return
    213214elseif (strcmpi(name,'ControlParameter')), enum=ControlParameterEnum(); return

  • issm/trunk/src/m/enum/TolXEnum.m

    r4043 r4047  
    99%      macro=TolXEnum()
    1010
    11 macro=244;
     11macro=245;

  • issm/trunk/src/m/enum/VariableDescriptorsEnum.m

    r4043 r4047  
    99%      macro=VariableDescriptorsEnum()
    1010
    11 macro=245;
     11macro=246;

  • issm/trunk/src/m/enum/VerboseEnum.m

    r4043 r4047  
    99%      macro=VerboseEnum()
    1010
    11 macro=246;
     11macro=247;

  • issm/trunk/src/m/enum/WaitOnLockEnum.m

    r4043 r4047  
    99%      macro=WaitOnLockEnum()
    1010
    11 macro=247;
     11macro=248;

  • issm/trunk/src/m/enum/YtsEnum.m

    r4043 r4047  
    99%      macro=YtsEnum()
    1010
    11 macro=248;
     11macro=249;
Note: See TracChangeset for help on using the changeset viewer.