Changeset 13283


Ignore:
Timestamp:
09/06/12 16:17:49 (13 years ago)
Author:
Eric.Larour
Message:

CHG: finished raw implementation to call ADOLC drivers.
The drivers are called in AutodiffDriversx.
We added parameters in Autodiff/CreateParametersAutodiff, such as the state variable xp.
Added a bunch of enums accordingly.
Hooked up all of this in the transient_core solutoin for now.

Of course, also added new module RequestedDependentsx, to declare which dependents are going
to be our diagnostics.

Location:
issm/trunk-jpl/src
Files:
8 added
21 edited

Legend:

Unmodified
Added
Removed

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

    r13268 r13283  
    1717        AutodiffIndependentsEnum,
    1818        AutodiffNumIndependentsEnum,
     19        AutodiffJacobianEnum,
     20        AutodiffXpEnum,
    1921        BalancethicknessSpcthicknessEnum,
    2022        BalancethicknessStabilizationEnum,

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

    r13268 r13283  
    329329                                        ./modules/RequestedOutputsx/RequestedOutputsx.h\
    330330                                        ./modules/RequestedOutputsx/RequestedOutputsx.cpp\
     331                                        ./modules/RequestedDependentsx/RequestedDependentsx.h\
     332                                        ./modules/RequestedDependentsx/RequestedDependentsx.cpp\
     333                                        ./modules/AutodiffDriversx/AutodiffDriversx.h\
     334                                        ./modules/AutodiffDriversx/AutodiffDriversx.cpp\
    331335                                        ./modules/ResetConstraintsx/ResetConstraintsx.h\
    332336                                        ./modules/ResetConstraintsx/ResetConstraintsx.cpp\

  • issm/trunk-jpl/src/c/classes/IoModel.cpp

    r13268 r13283  
    2626        this->fid=NULL;
    2727        this->data=NULL;
     28        this->independents=NULL;
    2829        this->constants=NULL;
    2930       
    … …  
    5556        this->data=xNew<IssmDouble*>(MaximumNumberOfEnums);
    5657        for(int i=0;i<MaximumNumberOfEnums;i++) this->data[i]=NULL;
     58
     59        /*Initialize array detecting whether data[i] is an independent AD mode variable: */
     60        this->independents=xNew<bool>(MaximumNumberOfEnums);
     61        for(int i=0;i<MaximumNumberOfEnums;i++) this->independents[i]=false;
    5762       
    5863        /*Initialize permanent data: */
    … …  
    8590
    8691        xDelete<IssmDouble*>(this->data);
     92        xDelete<bool>(this->independents);
    8793        xDelete<bool>(this->my_elements);
    8894        xDelete<bool>(this->my_nodes);
    … …  
    191197}
    192198/*}}}*/
    193 /*FUNCTION IoModel::DeleteData{{{*/
     199/*FUNCTION IoModel::DeleteData(int num,...){{{*/
    194200void  IoModel::DeleteData(int num,...){
    195201
    … …  
    205211                dataenum=va_arg(ap, int);
    206212                _assert_(dataenum<MaximumNumberOfEnums);
    207                 xDelete<IssmDouble>(this->data[dataenum]);
     213               
     214                /*do not erase independent variables for the AD mode computations!: */
     215                if (!this->independents[dataenum]) xDelete<IssmDouble>(this->data[dataenum]);
    208216        }
    209217        va_end(ap);
     218} /*}}}*/
     219/*FUNCTION IoModel::DeleteData(IssmDouble* {{{*/
     220void  IoModel::DeleteData(IssmDouble* vector, int dataenum){
     221
     222        /*do not erase independent variables for the AD mode computations!: */
     223        if(vector)if (!this->independents[dataenum]) xDelete<IssmDouble>(vector);
     224
    210225} /*}}}*/
    211226/*FUNCTION IoModel::FetchConstants{{{*/
    … …  
    508523        IssmPDouble   scalar;
    509524        int      code;
    510        
     525
    511526        /*Set file pointer to beginning of the data: */
    512527        fid=this->SetFilePointerToData(&code,NULL,data_enum);
    … …  
    538553        int   string_size;
    539554        int code=0;
    540        
     555
    541556        /*Set file pointer to beginning of the data: */
    542557        fid=this->SetFilePointerToData(&code,NULL,data_enum);
    … …  
    591606        int code=0;
    592607        int vector_type=0;
    593        
    594        
     608
     609
    595610        /*Set file pointer to beginning of the data: */
    596611        fid=this->SetFilePointerToData(&code,&vector_type,data_enum);
    … …  
    664679        int code=0;
    665680        int vector_type=0;
    666        
     681
    667682        /*Set file pointer to beginning of the data: */
    668683        fid=this->SetFilePointerToData(&code,&vector_type,data_enum);
    … …  
    698713                MPI_Bcast(matrix,M*N,MPI_DOUBLE,0,MPI_COMM_WORLD);
    699714                #endif
    700                 *pmatrix=xNew<IssmDouble>(M*N);
    701                 for (int i=0;i<M*N;++i) (*pmatrix)[i]=matrix[i];
     715
     716                if (this->independents[data_enum]){
     717                        /*this data has already been checked out! So cancel all that we've done here, and return
     718                         * the data[data_enum] directly: */
     719                        *pmatrix=this->data[data_enum];
     720                }
     721                else{
     722                        *pmatrix=xNew<IssmDouble>(M*N);
     723                        for (int i=0;i<M*N;++i) (*pmatrix)[i]=matrix[i];
     724                }
    702725                xDelete<IssmPDouble>(matrix);
    703726        }
    … …  
    936959               
    937960                dataenum=va_arg(ap, int);
     961               
     962                if (this->independents[dataenum]){
     963                        /*this data has already been checked out! Continue: */
     964                        continue;
     965                }
    938966
    939967                /*Some checks in debugging mode*/
    … …  
    11461174                        /*}}}*/
    11471175        }
    1148         /*Free ressources:*/
    1149         xDelete<IssmDouble>(IssmDoublevector);
     1176        /*Free ressources. Pay attention to not freeing an AD mode independent variable though!:*/
     1177        if (!this->independents[vector_enum] && !this->independents[default_vector_enum]) xDelete<IssmDouble>(IssmDoublevector);
    11501178        xDelete<char>(string);
    11511179}
     1180/*}}}*/
    11521181/*FUNCTION IoModel::LastIndex{{{*/
    11531182void IoModel::LastIndex(int *pindex){
    … …  
    11851214        *pindex=lastindex;
    11861215}
     1216/*}}}*/
    11871217/*FUNCTION IoModel::SetFilePointerToData{{{*/
    11881218FILE* IoModel::SetFilePointerToData(int* pcode,int* pvector_type, int data_enum){
    … …  
    12501280}
    12511281/*}}}*/
     1282/*FUNCTION IoModel::DeclareIndependents{{{*/
     1283void IoModel::DeclareIndependents(void){
     1284
     1285        bool autodiff=false;
     1286        int  dummy;
     1287        int i;
     1288        int  num_independents;
     1289        int* independents=NULL;
     1290
     1291        #ifdef _HAVE_ADOLC_
     1292        /*recover independent enums: */
     1293        this->Constant(&num_independents,AutodiffNumIndependentsEnum);
     1294        if(num_independents){
     1295                this->FetchData(&independents,&dummy,&dummy,AutodiffIndependentsEnum);
     1296
     1297                /*now go fetch the independent variables for each independent enum: */
     1298                for(i=0;i<num_independents;i++){
     1299                        this->FetchIndependent(independents[i]);
     1300                }
     1301                xDelete<int>(independents);
     1302        }
     1303        #else
     1304        /*if we asked for AD computations, we have a problem!: */
     1305        this->Constant(&autodiff,AutodiffIsautodiffEnum);
     1306        if(autodiff)_error_("Cannot carry out AD mode computations without support of ADOLC compiled in!");
     1307        #endif
     1308
     1309}
     1310/*}}}*/
     1311/*FUNCTION IoModel::FetchIndependent{{{*/
     1312void IoModel::FetchIndependent(int independent_enum){
     1313
     1314        #ifdef _HAVE_ADOLC_ //cannot come here unless you are running AD mode, from DeclaredIndependents:
     1315        extern int my_rank;
     1316        extern int num_procs;
     1317
     1318        /*output: */
     1319        int M,N;
     1320        IssmPDouble* buffer=NULL; //a buffer to read the data from disk
     1321        IssmDouble* matrix=NULL; //our independent variable
     1322        int code=0;
     1323        int vector_type=0;
     1324       
     1325        /*Set file pointer to beginning of the data: */
     1326        fid=this->SetFilePointerToData(&code,&vector_type,independent_enum);
     1327        if((code!=5) && (code!=6) && (code!=7))_error_("expecting a IssmDouble, integer or boolean matrix for enum " << EnumToStringx(independent_enum));
     1328       
     1329        /*Now fetch: */
     1330
     1331        /*We have to read a matrix from disk. First read the dimensions of the matrix, then the whole matrix: */
     1332        /*numberofelements: */
     1333        if(my_rank==0){ 
     1334                if(fread(&M,sizeof(int),1,fid)!=1) _error_("could not read number of rows for matrix ");
     1335        }
     1336        #ifdef _HAVE_MPI_
     1337        MPI_Bcast(&M,1,MPI_INT,0,MPI_COMM_WORLD);
     1338        #endif
     1339
     1340        if(my_rank==0){ 
     1341                if(fread(&N,sizeof(int),1,fid)!=1) _error_("could not read number of columns for matrix ");
     1342        }
     1343        #ifdef _HAVE_MPI_
     1344        MPI_Bcast(&N,1,MPI_INT,0,MPI_COMM_WORLD);
     1345        #endif
     1346
     1347        /*Now allocate matrix: */
     1348        if(M*N){
     1349                buffer=xNew<IssmPDouble>(M*N);
     1350                matrix=xNew<IssmDouble>(M*N);
     1351
     1352                /*Read matrix on node 0, then broadcast: */
     1353                if(my_rank==0){ 
     1354                        if(fread(buffer,M*N*sizeof(IssmPDouble),1,fid)!=1) _error_("could not read matrix ");
     1355                       
     1356                        /*Now, before we even broadcast this to other nodes, declare the whole matrix as a independent variable!: */
     1357                        for (int i=0;i<M*N;++i) matrix[i]<<=buffer[i];  /*we use the <<= ADOLC overloaded operator to declare the independency*/
     1358                }
     1359                #ifdef _HAVE_MPI_
     1360                MPI_Bcast(matrix,M*N,MPI_DOUBLE,0,MPI_COMM_WORLD);
     1361                #endif
     1362               
     1363                xDelete<IssmPDouble>(buffer);
     1364        }
     1365        else _error_("cannot declare the independent variable " << EnumToStringx(independent_enum) <<  "if it's empty!");
     1366
     1367        /*Ok, we are almost done. Matrix is now a independent matrix. We don't want this matrix to be fetched again in the
     1368         *future, which would effectively write over the independency in the ADOLC tape! So we are going to keep track of this
     1369         independent matrix inthe iomodel->data[independent_enum] data slot: */
     1370        this->data[independent_enum]=matrix;
     1371        this->independents[independent_enum]=true;
     1372        #endif
     1373}
     1374/*}}}*/

  • issm/trunk-jpl/src/c/classes/IoModel.h

    r13268 r13283  
    1919
    2020        private:
     21                Parameters  *constants;   //this dataset holds all IssmDouble, int, bool and char *parameters read in from the input file.*
     22       
     23        public:
    2124                IssmDouble **data;        //this dataset holds temporary data, memory intensive.
    22                 Parameters  *constants;   //this dataset holds all IssmDouble, int, bool and char *parameters read in from the input file.*
    23 
    24         public:
     25               
    2526                /*This data needs to stay memory resident at all time, even if it's memory intensive: */
    2627                FILE *fid;         //pointer to input file
    … …  
    3132                int  *numbernodetoelementconnectivity;
    3233
     34
    3335                /*Data to synchronize through low level object drivers: */
    3436                int nodecounter;         //keep track of how many nodes are being created in each analysis type
    3537                int loadcounter;         //keep track of how many loads are being created in each analysis type
    3638                int constraintcounter;   //keep track of how many constraints are being created in each analysis type
     39               
     40                /*for AD mode: to keep track of our independent variables we fetch:*/
     41                bool* independents;
    3742
    3843                /*Methods*/
    … …  
    5055                IssmDouble *Data(int dataenum);
    5156                void        DeleteData(int num,...);
     57                void        DeleteData(IssmDouble* vector, int dataenum);
    5258                void        FetchConstants(void);
    5359                void        FetchData(bool* pboolean,int data_enum);
    … …  
    6470                void        LastIndex(int *pindex);
    6571                FILE*       SetFilePointerToData(int* pcode,int* pvector_type, int data_enum);
    66                 #ifdef _HAVE_ADOLC_
    67                 void        FetchIndependentVariable(IssmDouble**  pscalarmatrix,int* pM,int* pN,int data_enum);
    68                 #endif
     72                void        DeclareIndependents(void);
     73                void        FetchIndependent(int dependent_enum);
    6974};
    7075

  • issm/trunk-jpl/src/c/modules/EnumToStringx/EnumToStringx.cpp

    r13268 r13283  
    2222                case AutodiffIndependentsEnum : return "AutodiffIndependents";
    2323                case AutodiffNumIndependentsEnum : return "AutodiffNumIndependents";
     24                case AutodiffJacobianEnum : return "AutodiffJacobian";
     25                case AutodiffXpEnum : return "AutodiffXp";
    2426                case BalancethicknessSpcthicknessEnum : return "BalancethicknessSpcthickness";
    2527                case BalancethicknessStabilizationEnum : return "BalancethicknessStabilization";

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/Autodiff/CreateParametersAutodiff.cpp

    r13268 r13283  
    1414void CreateParametersAutodiff(Parameters** pparameters,IoModel* iomodel,int solution_type,int analysis_type){
    1515       
    16         int         i;
     16        int         i,j;
    1717        Parameters *parameters       = NULL;
    1818        bool        autodiff_analysis;
    … …  
    2121        int*        independents       = NULL;
    2222        int         num_independents;
     23        int         numberofvertices;
     24        IssmDouble* xp=NULL;
    2325       
    2426        /*Get parameters: */
    … …  
    2729        /*retrieve some parameters: */
    2830        iomodel->Constant(&autodiff_analysis,AutodiffIsautodiffEnum);
     31       
     32        if(autodiff_analysis){
    2933
    30         if(autodiff_analysis){
     34                iomodel->Constant(&num_independents,AutodiffNumIndependentsEnum);
     35                iomodel->Constant(&num_dependents,AutodiffNumDependentsEnum);
     36                iomodel->Constant(&numberofvertices,MeshNumberofverticesEnum);
    3137
    3238                /*recover dependents: */
    3339                parameters->AddObject(iomodel->CopyConstantObject(AutodiffNumDependentsEnum));
    34                 iomodel->FetchData(&dependents,NULL,&num_dependents,AutodiffDependentsEnum);
    35                 parameters->AddObject(new IntVecParam(AutodiffDependentsEnum,dependents,num_dependents));
    36        
     40                if(num_dependents){
     41                        iomodel->FetchData(&dependents,NULL,&num_dependents,AutodiffDependentsEnum);
     42                        parameters->AddObject(new IntVecParam(AutodiffDependentsEnum,dependents,num_dependents));
     43                }
     44
    3745                /*recover independents: */
    3846                parameters->AddObject(iomodel->CopyConstantObject(AutodiffNumIndependentsEnum));
    39                 iomodel->FetchData(&independents,NULL,&num_independents,AutodiffIndependentsEnum);
    40                 parameters->AddObject(new IntVecParam(AutodiffIndependentsEnum,independents,num_independents));
     47                if(num_independents){
     48                        iomodel->FetchData(&independents,NULL,&num_independents,AutodiffIndependentsEnum);
     49                        parameters->AddObject(new IntVecParam(AutodiffIndependentsEnum,independents,num_independents));
    4150
     51                        /*Build state vector, value at which we compute our gradients of dependent variables in adolc: the xp vector  */
     52                        xp=xNew<IssmDouble>(num_independents*numberofvertices);
     53                        for(i=0;i<num_independents;i++){
     54                                IssmDouble* values=iomodel->data[independents[i]];
     55                                for(j=0;j<numberofvertices;j++){
     56                                        xp[i*numberofvertices+j]=values[j];
     57                                }
     58                        }
     59                        parameters->AddObject(new DoubleVecParam(AutodiffXpEnum,xp,num_independents*numberofvertices));
     60                }
     61
     62                /*Assign output pointer: */
     63                *pparameters=parameters;
    4264        }
    43 
    44         /*Assign output pointer: */
    45         *pparameters=parameters;
    4665}

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/Control/CreateParametersControl.cpp

    r13073 r13283  
    6464
    6565                xDelete<int>(control_type);
    66                 xDelete<IssmDouble>(cm_responses);
    67                 xDelete<IssmDouble>(cm_jump);
    68                 xDelete<IssmDouble>(optscal);
    69                 xDelete<IssmDouble>(maxiter);
     66                iomodel->DeleteData(cm_responses,InversionCostFunctionsEnum);
     67                iomodel->DeleteData(cm_jump,InversionStepThresholdEnum);
     68                iomodel->DeleteData(optscal,InversionGradientScalingEnum);
     69                iomodel->DeleteData(maxiter,InversionMaxiterPerStepEnum);
    7070        }
    7171

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateParameters.cpp

    r13277 r13283  
    9292        parameters->AddObject(iomodel->CopyConstantObject(SurfaceforcingsIsdelta18oEnum));
    9393        parameters->AddObject(iomodel->CopyConstantObject(SurfaceforcingsIssmbgradientsEnum));
    94         parameters->AddObject(iomodel->CopyConstantObject(AutodiffIsautodiffEnum));
    9594
    9695        iomodel->Constant(&isdelta18o,SurfaceforcingsIsdelta18oEnum);
    … …  
    105104                for(i=0;i<M;i++) temp[M+i]=yts*temp[M+i];
    106105                parameters->AddObject(new TransientParam(SurfaceforcingsDelta18oEnum,&temp[0],&temp[M],M));
    107                 xDelete<IssmDouble>(temp);
     106                iomodel->DeleteData(temp,SurfaceforcingsDelta18oEnum);
    108107
    109108                iomodel->FetchData(&temp,&N,&M,SurfaceforcingsDelta18oSurfaceEnum); _assert_(N==2);
    110109                for(i=0;i<M;i++) temp[M+i]=yts*temp[M+i];
    111110                parameters->AddObject(new TransientParam(SurfaceforcingsDelta18oSurfaceEnum,&temp[0],&temp[M],M));
    112                 xDelete<IssmDouble>(temp);
     111                iomodel->DeleteData(temp,SurfaceforcingsDelta18oSurfaceEnum);
    113112        }
    114113
    … …  
    127126        parameters->AddObject(new IntParam(DiagnosticNumRequestedOutputsEnum,numoutputs));
    128127        if(numoutputs)parameters->AddObject(new IntVecParam(DiagnosticRequestedOutputsEnum,requestedoutputs,numoutputs));
    129         xDelete<IssmDouble>(requestedoutputs);
     128        iomodel->DeleteData(requestedoutputs,DiagnosticRequestedOutputsEnum);
    130129        iomodel->FetchData(&requestedoutputs,&numoutputs,NULL,TransientRequestedOutputsEnum);
    131130        parameters->AddObject(new IntParam(TransientNumRequestedOutputsEnum,numoutputs));
    132131        if(numoutputs)parameters->AddObject(new IntVecParam(TransientRequestedOutputsEnum,requestedoutputs,numoutputs));
    133         xDelete<IssmDouble>(requestedoutputs);
     132        iomodel->DeleteData(requestedoutputs,TransientRequestedOutputsEnum);
    134133        iomodel->FetchData(&requestedoutputs,&numoutputs,NULL,SteadystateRequestedOutputsEnum);
    135134        parameters->AddObject(new IntParam(SteadystateNumRequestedOutputsEnum,numoutputs));
    136135        if(numoutputs)parameters->AddObject(new IntVecParam(SteadystateRequestedOutputsEnum,requestedoutputs,numoutputs));
    137         xDelete<IssmDouble>(requestedoutputs);
     136        iomodel->DeleteData(requestedoutputs,SteadystateRequestedOutputsEnum);
    138137       
    139138        /*Before returning, create parameters in case we are running Qmu or control types runs: */

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateConstraintsDiagnosticHoriz.cpp

    r13073 r13283  
    326326         
    327327        /*Free data: */
    328         xDelete<IssmDouble>(spcvx);
    329         xDelete<IssmDouble>(spcvy);
    330         xDelete<IssmDouble>(spcvz);
    331         xDelete<IssmDouble>(nodeonmacayeal);
    332         xDelete<IssmDouble>(nodeonpattyn);
    333         xDelete<IssmDouble>(nodeonstokes);
    334         xDelete<IssmDouble>(nodeonicesheet);
    335         xDelete<IssmDouble>(nodeonbed);
    336         xDelete<IssmDouble>(vertices_type);
    337         xDelete<IssmDouble>(surface);
    338         xDelete<IssmDouble>(z);
     328        iomodel->DeleteData(spcvx,DiagnosticSpcvxEnum);
     329        iomodel->DeleteData(spcvy,DiagnosticSpcvyEnum);
     330        iomodel->DeleteData(spcvz,DiagnosticSpcvzEnum);
     331        iomodel->DeleteData(nodeonmacayeal,FlowequationBordermacayealEnum);
     332        if(dim==3)iomodel->DeleteData(nodeonpattyn,FlowequationBorderpattynEnum);
     333        if(dim==3)iomodel->DeleteData(nodeonstokes,FlowequationBorderstokesEnum);
     334        if(dim==3)iomodel->DeleteData(nodeonbed,MeshVertexonbedEnum);
     335        if(dim==3)iomodel->DeleteData(nodeonicesheet,MaskVertexongroundediceEnum);
     336        iomodel->DeleteData(vertices_type,FlowequationVertexEquationEnum);
     337        iomodel->DeleteData(surface,SurfaceEnum);
     338        iomodel->DeleteData(z,MeshZEnum);
    339339
    340340        /*Free resources:*/

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/DiagnosticHoriz/CreateLoadsDiagnosticHoriz.cpp

    r13073 r13283  
    129129        /*Free data: */
    130130        iomodel->DeleteData(3,DiagnosticIcefrontEnum,ThicknessEnum,BedEnum);
    131         xDelete<IssmDouble>(elements_type);
    132         xDelete<IssmDouble>(pressureload);
     131        iomodel->DeleteData(elements_type,FlowequationElementEquationEnum);
     132        iomodel->DeleteData(pressureload,DiagnosticIcefrontEnum);
     133
    133134
    134135        /*Create Penpair for penalties: */
    … …  
    153154
    154155        /*free ressources: */
    155         xDelete<IssmDouble>(penalties);
     156        iomodel->DeleteData(penalties,DiagnosticVertexPairingEnum);
    156157
    157158        /*Create Riffront loads for rifts: */

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/ElementsAndVerticesPartitioning.cpp

    r12832 r13283  
    7676
    7777        /*Free elements and elements2d: */
    78         xDelete<IssmDouble>(elements);
    79         xDelete<IssmDouble>(elements2d);
     78        iomodel->DeleteData(elements,MeshElementsEnum);
     79        iomodel->DeleteData(elements2d,MeshElements2dEnum);
    8080
    8181        /*Deal with rifts, they have to be included into one partition only, not several: */
    … …  
    8787                        epart[el2]=epart[el1]; //ensures that this pair of elements will be in the same partition, as well as the corresponding vertices;
    8888                }
    89                 xDelete<IssmDouble>(riftinfo);
     89                iomodel->DeleteData(riftinfo,RiftsRiftstructEnum);
    9090        }
    9191
    … …  
    119119        }//for (i=0;i<numberofelements;i++)
    120120        /*Free data : */
    121         xDelete<IssmDouble>(elements);
     121        iomodel->DeleteData(elements,MeshElementsEnum);
    122122
    123123        /*We might have vertex_pairing in which case, some vertices have to be cloned:
    … …  
    130130                }
    131131        }
    132         xDelete<IssmDouble>(vertex_pairing);
     132        iomodel->DeleteData(vertex_pairing,DiagnosticVertexPairingEnum);
    133133        iomodel->FetchData(&vertex_pairing,&numvertex_pairing,NULL,PrognosticVertexPairingEnum);
    134134        for(i=0;i<numvertex_pairing;i++){
    … …  
    137137                }
    138138        }
    139         xDelete<IssmDouble>(vertex_pairing);
     139        iomodel->DeleteData(vertex_pairing,PrognosticVertexPairingEnum);
    140140
    141141        /*Free ressources:*/

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/Enthalpy/CreateConstraintsEnthalpy.cpp

    r13073 r13283  
    110110
    111111        /*Free ressources:*/
    112         xDelete<IssmDouble>(spcvector);
     112        iomodel->DeleteData(spcvector,ThermalSpctemperatureEnum);
    113113        xDelete<IssmDouble>(times);
    114114        xDelete<IssmDouble>(values);

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/ModelProcessorx.cpp

    r13277 r13283  
    2020
    2121        int   i,analysis_type,dim,verbose;
    22         bool  isthermal,isprognostic,isdiagnostic,isgroundingline,isenthalpy;
     22        bool  isthermal,isprognostic,isdiagnostic,isgroundingline,isenthalpy,autodiff;
    2323       
    2424        /*output: */
    … …  
    4343        iomodel->Constant(&isdiagnostic,TransientIsdiagnosticEnum);
    4444        iomodel->Constant(&isgroundingline,TransientIsgroundinglineEnum);
    45        
     45        iomodel->Constant(&autodiff,AutodiffIsautodiffEnum);
     46
     47        /*If we are running in AD mode, we need to declare our independent variables now, before
     48         *and prevent them from being erased during successive calls to iomodel->FetchData and
     49         iomodel->DeleteData:*/
     50        if(autodiff)iomodel->DeclareIndependents();
     51
    4652        SetVerbosityLevel(verbose);
    4753

  • issm/trunk-jpl/src/c/modules/ModelProcessorx/Prognostic/CreateLoadsPrognostic.cpp

    r12832 r13283  
    9090
    9191        /*free ressources: */
    92         xDelete<IssmDouble>(vertex_pairing);
    93         xDelete<IssmDouble>(nodeonbed);
     92        iomodel->DeleteData(vertex_pairing,PrognosticVertexPairingEnum);
     93        iomodel->DeleteData(nodeonbed,MeshVertexonbedEnum);
    9494
    9595        /*Assign output pointer: */

  • issm/trunk-jpl/src/c/modules/StringToEnumx/StringToEnumx.cpp

    r13268 r13283  
    2323              else if (strcmp(name,"AutodiffIndependents")==0) return AutodiffIndependentsEnum;
    2424              else if (strcmp(name,"AutodiffNumIndependents")==0) return AutodiffNumIndependentsEnum;
     25              else if (strcmp(name,"AutodiffJacobian")==0) return AutodiffJacobianEnum;
     26              else if (strcmp(name,"AutodiffXp")==0) return AutodiffXpEnum;
    2527              else if (strcmp(name,"BalancethicknessSpcthickness")==0) return BalancethicknessSpcthicknessEnum;
    2628              else if (strcmp(name,"BalancethicknessStabilization")==0) return BalancethicknessStabilizationEnum;
    … …  
    136138              else if (strcmp(name,"MeshVertexonbed")==0) return MeshVertexonbedEnum;
    137139              else if (strcmp(name,"MeshVertexonsurface")==0) return MeshVertexonsurfaceEnum;
    138               else if (strcmp(name,"MeshX")==0) return MeshXEnum;
    139               else if (strcmp(name,"MeshY")==0) return MeshYEnum;
    140140         else stage=2;
    141141   }
    142142   if(stage==2){
    143               if (strcmp(name,"MeshZ")==0) return MeshZEnum;
     143              if (strcmp(name,"MeshX")==0) return MeshXEnum;
     144              else if (strcmp(name,"MeshY")==0) return MeshYEnum;
     145              else if (strcmp(name,"MeshZ")==0) return MeshZEnum;
    144146              else if (strcmp(name,"MiscellaneousName")==0) return MiscellaneousNameEnum;
    145147              else if (strcmp(name,"PrognosticHydrostaticAdjustment")==0) return PrognosticHydrostaticAdjustmentEnum;
    … …  
    259261              else if (strcmp(name,"Materials")==0) return MaterialsEnum;
    260262              else if (strcmp(name,"Nodes")==0) return NodesEnum;
    261               else if (strcmp(name,"Parameters")==0) return ParametersEnum;
    262               else if (strcmp(name,"Vertices")==0) return VerticesEnum;
    263263         else stage=3;
    264264   }
    265265   if(stage==3){
    266               if (strcmp(name,"Results")==0) return ResultsEnum;
     266              if (strcmp(name,"Parameters")==0) return ParametersEnum;
     267              else if (strcmp(name,"Vertices")==0) return VerticesEnum;
     268              else if (strcmp(name,"Results")==0) return ResultsEnum;
    267269              else if (strcmp(name,"AdolcParam")==0) return AdolcParamEnum;
    268270              else if (strcmp(name,"BoolInput")==0) return BoolInputEnum;
    … …  
    382384              else if (strcmp(name,"VzMacAyeal")==0) return VzMacAyealEnum;
    383385              else if (strcmp(name,"VzPattyn")==0) return VzPattynEnum;
    384               else if (strcmp(name,"VzPicard")==0) return VzPicardEnum;
    385               else if (strcmp(name,"VzStokes")==0) return VzStokesEnum;
    386386         else stage=4;
    387387   }
    388388   if(stage==4){
    389               if (strcmp(name,"VxMesh")==0) return VxMeshEnum;
     389              if (strcmp(name,"VzPicard")==0) return VzPicardEnum;
     390              else if (strcmp(name,"VzStokes")==0) return VzStokesEnum;
     391              else if (strcmp(name,"VxMesh")==0) return VxMeshEnum;
    390392              else if (strcmp(name,"VyMesh")==0) return VyMeshEnum;
    391393              else if (strcmp(name,"VzMesh")==0) return VzMeshEnum;

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

    r12927 r13283  
    88/*Modules: */
    99#include "./AddExternalResultx/AddExternalResultx.h"
     10#include "./AutodiffDriversx/AutodiffDriversx.h"
    1011#include "./AverageFilterx/AverageFilterx.h"
    1112#include "./AverageOntoPartitionx/AverageOntoPartitionx.h"
    … …  
    103104#include "./Reducevectorgtofx/Reducevectorgtofx.h"
    104105#include "./RequestedOutputsx/RequestedOutputsx.h"
     106#include "./RequestedDependentsx/RequestedDependentsx.h"
    105107#include "./ResetConstraintsx/ResetConstraintsx.h"
    106108#include "./ResetCoordinateSystemx/ResetCoordinateSystemx.h"

  • issm/trunk-jpl/src/c/solutions/issm.cpp

    r13277 r13283  
    8484        /*out of solution_type, figure out types of analyses needed in the femmodel: */
    8585        AnalysisConfiguration(&analyses,&numanalyses,solution_type);
    86        
     86
     87        /*before we create the model, start the trace on for AD mode: */
     88        #ifdef _HAVE_ADOLC_
     89        trace_on(1);
     90        #endif
     91
    8792        /*Create femmodel, using input file: */
    8893        #ifdef _HAVE_MPI_
    … …  
    124129        theAdolcEDF_p->GetParameterValue().myEDF_for_solverx_p=reg_ext_fct(EDF_for_solverx);
    125130        femmodel->parameters->AddObject(theAdolcEDF_p);
    126         if(autodiff) trace_on(1);
    127131        #else
    128132        if(autodiff) _error_("ISSM was not compiled with ADOLC support, cannot carry out autodiff analysis!");
    … …  
    166170        }
    167171
     172        if(autodiff){
     173                trace_off();
     174                AutodiffDriversx(femmodel->elements, femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,femmodel->results);
     175        }
     176
    168177
    169178        #ifdef _HAVE_MPI_
    … …  
    172181        finish_core=(IssmPDouble)clock();
    173182        #endif
    174 
     183       
    175184        _pprintLine_("write results to disk:");
    176185        OutputResultsx(femmodel->elements, femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,femmodel->results);
    … …  
    179188        #ifdef _HAVE_ADOLC_
    180189        if(autodiff){
    181                 trace_off();
    182190                tapestats(1,tape_stats); //reading of tape statistics
    183191                _pprintLine_("   ADOLC statistics: ");

  • issm/trunk-jpl/src/c/solutions/transient_core.cpp

    r13056 r13283  
    152152        }
    153153
     154        RequestedDependentsx(femmodel->results,femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters);
     155       
    154156        /*Free ressources:*/
    155157        xDelete<int>(requested_outputs);

  • issm/trunk-jpl/src/m/classes/autodiff.m

    r13271 r13283  
    77        properties (SetAccess=public)
    88                isautodiff = false;
    9                 dependents = {''};
    10                 independents = {''};
     9                dependents = {};
     10                independents = {};
    1111        end
    1212        methods
    … …  
    2323                end % }}}
    2424                function md = checkconsistency(obj,md,solution,analyses) % {{{
     25
     26                %Early return
     27                if ~obj.isautodiff, return; end
     28
    2529
    2630                end % }}}
    … …  
    4044                        %process dependent variables
    4145                        num_dependents=numel(obj.dependents);
    42                         data=zeros(1,num_dependents);
    43                         for i=1:num_dependents,
    44                                 data(i)=StringToEnum(obj.dependents{i});
     46                        WriteData(fid,'data',num_dependents,'enum',AutodiffNumDependentsEnum(),'format','Integer');
     47                       
     48                        if(num_dependents),
     49                                data=zeros(1,num_dependents);
     50                                for i=1:num_dependents,
     51                                        data(i)=StringToEnum(obj.dependents{i});
     52                                end
     53                                WriteData(fid,'data',data,'enum',AutodiffDependentsEnum(),'format','DoubleMat','mattype',3);
    4554                        end
    46                         WriteData(fid,'data',data,'enum',AutodiffDependentsEnum(),'format','DoubleMat','mattype',3);
    47                         WriteData(fid,'data',num_dependents,'enum',AutodiffNumDependentsEnum(),'format','Integer');
    4855                       
    4956                        %process independent variables
    5057                        num_independents=numel(obj.independents);
    51                         data=zeros(1,num_independents);
    52                         for i=1:num_independents,
    53                                 data(i)=StringToEnum(obj.independents{i});
     58                        WriteData(fid,'data',num_independents,'enum',AutodiffNumIndependentsEnum(),'format','Integer');
     59                       
     60                        if(num_independents)
     61                                data=zeros(1,num_independents);
     62                                for i=1:num_independents,
     63                                        data(i)=StringToEnum(obj.independents{i});
     64                                end
     65                                WriteData(fid,'data',data,'enum',AutodiffIndependentsEnum(),'format','DoubleMat','mattype',3);
    5466                        end
    55                         WriteData(fid,'data',data,'enum',AutodiffIndependentsEnum(),'format','DoubleMat','mattype',3);
    56                         WriteData(fid,'data',num_independents,'enum',AutodiffNumIndependentsEnum(),'format','Integer');
    5767
    5868                end % }}}

  • issm/trunk-jpl/src/m/enum/EnumDefinitions.py

    r13268 r13283  
    5959        return StringToEnum('AutodiffNumIndependents')[0]
    6060
     61def AutodiffJacobianEnum():
     62        """
     63        AUTODIFFJACOBIANENUM - Enum of AutodiffJacobian
     64
     65           Usage:
     66              macro=AutodiffJacobianEnum()
     67        """
     68
     69        return StringToEnum('AutodiffJacobian')[0]
     70
     71def AutodiffXpEnum():
     72        """
     73        AUTODIFFXPENUM - Enum of AutodiffXp
     74
     75           Usage:
     76              macro=AutodiffXpEnum()
     77        """
     78
     79        return StringToEnum('AutodiffXp')[0]
     80
    6181def BalancethicknessSpcthicknessEnum():
    6282        """
    … …  
    47174737        """
    47184738
    4719         return 470
    4720 
     4739        return 472
     4740

  • issm/trunk-jpl/src/m/enum/MaximumNumberOfEnums.m

    r13268 r13283  
    99%      macro=MaximumNumberOfEnums()
    1010
    11 macro=470;
     11macro=472;
Note: See TracChangeset for help on using the changeset viewer.