Changeset 3840


Ignore:
Timestamp:
05/19/10 09:57:14 (15 years ago)
Author:
Mathieu Morlighem
Message:

some improvements and bug fix

Location:
issm/trunk/src/c/objects
Files:
6 edited

Legend:

Unmodified
Added
Removed

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

    r3833 r3840  
    545545}
    546546/*}}}*/
    547 /*FUNCTION Penta::GetSolutionFromInputs(Vec solution,  int analysis_type,int sub_analysis_type);{{{1*/
    548 void  Penta::GetSolutionFromInputs(Vec solution,  int analysis_type,int sub_analysis_type){
    549         ISSMERROR(" not supported yet!");
     547/*FUNCTION Penta::GetSolutionFromInputs(Vec solution,int analysis_type,int sub_analysis_type){{{1*/
     548void  Penta::GetSolutionFromInputs(Vec solution,int analysis_type,int sub_analysis_type){
     549        /*Just branch to the correct UpdateInputsFromSolution generator, according to the type of analysis we are carrying out: */
     550        if (analysis_type==DiagnosticAnalysisEnum){
     551                if (sub_analysis_type==HorizAnalysisEnum){
     552                        GetSolutionFromInputsDiagnosticHoriz(solution,analysis_type,sub_analysis_type);
     553                }
     554                else ISSMERROR("%s%i%s\n","sub_analysis: ",sub_analysis_type," not supported yet");
     555        }
     556        else{
     557                ISSMERROR("%s%i%s\n","analysis: ",analysis_type," not supported yet");
     558        }
     559}
     560/*}}}*/
     561/*FUNCTION Penta::GetSolutionFromInputsDiagnosticHoriz(Vec solution,int analysis_type,int sub_analysis_type){{{1*/
     562void  Penta::GetSolutionFromInputsDiagnosticHoriz(Vec solution,int analysis_type,int sub_analysis_type){
     563
     564        int i;
     565
     566        const int    numvertices=6;
     567        const int    numdofpervertex=2;
     568        const int    numdof=numdofpervertex*numvertices;
     569        double       gauss[numvertices][4]={{1,0,0,-1},{0,1,0,-1},{0,0,1,-1},{1,0,0,1},{0,1,0,1},{0,0,1,1}};
     570
     571        int          doflist[numdof];
     572        double       values[numdof];
     573        double       vx;
     574        double       vy;
     575
     576        int          dummy;
     577
     578        /*Get dof list: */
     579        GetDofList(&doflist[0],&dummy);
     580
     581        /*Ok, we have vx and vy in values, fill in vx and vy arrays: */
     582        /*P1 element only for now*/
     583        for(i=0;i<numvertices;i++){
     584
     585                /*Recover vx and vy*/
     586                inputs->GetParameterValue(&vx,&gauss[i][0],VxEnum);
     587                inputs->GetParameterValue(&vy,&gauss[i][0],VyEnum);
     588                values[i*numdofpervertex+0]=vx;
     589                values[i*numdofpervertex+1]=vy;
     590        }
     591
     592        /*Add value to global vector*/
     593        VecSetValues(solution,numdof,doflist,(const double*)values,INSERT_VALUES);
     594
    550595}
    551596/*}}}*/
    … …  
    11691214        double* area_gauss_weights  =  NULL;
    11701215        double* vert_gauss_weights  =  NULL;
    1171         double  gaussgrids[numgrids][numgrids]={{1,0,0,0},{0,1,0,0},{0,0,1,0},{1,0,0,1},{0,1,0,1},{0,0,1,1}};
     1216        double  gaussgrids[numgrids][4]={{1,0,0,-1},{0,1,0,-1},{0,0,1,-1},{1,0,0,1},{0,1,0,1},{0,0,1,1}};
    11721217
    11731218        /* specific gaussian point: */
    … …  
    30643109
    30653110        const int numgrids=6;
    3066         double  gaussgrids[numgrids][numgrids]={{1,0,0,0},{0,1,0,0},{0,0,1,0},{1,0,0,1},{0,1,0,1},{0,0,1,1}};
     3111        double  gaussgrids[numgrids][4]={{1,0,0,-1},{0,1,0,-1},{0,0,1,-1},{1,0,0,1},{0,1,0,1},{0,0,1,1}};
    30673112       
    30683113        inputs->GetParameterValues(bedlist,&gaussgrids[0][0],6,BedEnum);
    … …  
    39894034
    39904035        const int numgrids=6;
    3991         double  gaussgrids[numgrids][numgrids]={{1,0,0,0},{0,1,0,0},{0,0,1,0},{1,0,0,1},{0,1,0,1},{0,0,1,1}};
     4036        double  gaussgrids[numgrids][4]={{1,0,0,-1},{0,1,0,-1},{0,0,1,-1},{1,0,0,1},{0,1,0,1},{0,0,1,1}};
    39924037        inputs->GetParameterValues(thickness_list,&gaussgrids[0][0],6,ThicknessEnum);
    39934038

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

    r3821 r3840  
    6363                void  CreatePVector(Vec pg,  int analysis_type,int sub_analysis_type);
    6464                void  GetSolutionFromInputs(Vec solution,  int analysis_type,int sub_analysis_type);
     65                void  GetSolutionFromInputsDiagnosticHoriz(Vec solution,int analysis_type,int sub_analysis_type);
    6566                void  GetDofList(int* doflist,int* pnumberofdofs);
    6667                void  GetDofList1(int* doflist);

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

    r3834 r3840  
    483483}
    484484/*}}}*/
    485 /*FUNCTION Tria::GetSolutionFromInputs(Vec solution,  int analysis_type,int sub_analysis_type){{{1*/
    486 void  Tria::GetSolutionFromInputs(Vec solution,  int analysis_type,int sub_analysis_type){
     485/*FUNCTION Tria::GetSolutionFromInputs(Vec solution,int analysis_type,int sub_analysis_type){{{1*/
     486void  Tria::GetSolutionFromInputs(Vec solution,int analysis_type,int sub_analysis_type){
    487487        /*Just branch to the correct UpdateInputsFromSolution generator, according to the type of analysis we are carrying out: */
    488488        if (analysis_type==DiagnosticAnalysisEnum){
    … …  
    497497}
    498498/*}}}*/
    499 /*FUNCTION Tria::GetSolutionFromInputsDiagnosticHoriz(Vec solution, int analysis_type,int sub_analysis_type){{{1*/
     499/*FUNCTION Tria::GetSolutionFromInputsDiagnosticHoriz(Vec solution,int analysis_type,int sub_analysis_type){{{1*/
    500500void  Tria::GetSolutionFromInputsDiagnosticHoriz(Vec solution,int analysis_type,int sub_analysis_type){
    501501

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

    r3830 r3840  
    152152/*}}}*/
    153153/*FUNCTION PentaVertexInput::GetParameterValue(double* pvalue,double* gauss){{{1*/
    154 void PentaVertexInput::GetParameterValue(double* pvalue,double* gauss){ISSMERROR(" not supported yet!");}
     154void PentaVertexInput::GetParameterValue(double* pvalue,double* gauss){
     155        /*P1 interpolation on Gauss point*/
     156
     157        /*intermediary*/
     158        double l1l6[6];
     159
     160        /*nodal functions: */
     161        GetNodalFunctionsP1(&l1l6[0],gauss);
     162
     163        /*Assign output pointers:*/
     164        *pvalue=l1l6[0]*values[0]+l1l6[1]*values[1]+l1l6[2]*values[2]+l1l6[3]*values[3]+l1l6[4]*values[4]+l1l6[5]*values[5];
     165
     166}
    155167/*}}}*/
    156168/*FUNCTION PentaVertexInput::GetParameterValue(double* pvalue,double* gauss,double defaultvalue){{{1*/
    … …  
    158170/*}}}*/
    159171/*FUNCTION PentaVertexInput::GetParameterValues(double* values,double* gauss_pointers, int numgauss){{{1*/
    160 void PentaVertexInput::GetParameterValues(double* values,double* gauss_pointers, int numgauss){ISSMERROR(" not supported yet!");}
     172void PentaVertexInput::GetParameterValues(double* values,double* gauss_pointers, int numgauss){
     173        /*It is assumed that output values has been correctly allocated*/
     174
     175        int i,j;
     176        double gauss[4];
     177
     178        for (i=0;i<numgauss;i++){
     179
     180                /*Get current Gauss point coordinates*/
     181                for (j=0;j<4;j++) gauss[j]=gauss_pointers[i*4+j];
     182
     183                /*Assign parameter value*/
     184                GetParameterValue(&values[i],&gauss[0]);
     185        }
     186}
    161187/*}}}*/
    162188/*FUNCTION PentaVertexInput::GetParameterDerivativeValue(double* derivativevalues, double* xyz_list, double* gauss){{{1*/
    163 void PentaVertexInput::GetParameterDerivativeValue(double* derivativevalues, double* xyz_list, double* gauss){ISSMERROR(" not supported yet!");}
     189void PentaVertexInput::GetParameterDerivativeValue(double* p, double* xyz_list, double* gauss){
     190        /*From grid values of parameter p (p_list[0], p_list[1], p_list[2], p_list[3], p_list[4] and p_list[4]), return parameter derivative value at gaussian point specified by gauss_coord:
     191         *   dp/dx=p_list[0]*dh1/dx+p_list[1]*dh2/dx+p_list[2]*dh3/dx+p_list[3]*dh4/dx+p_list[4]*dh5/dx+p_list[5]*dh6/dx;
     192         *   dp/dy=p_list[0]*dh1/dy+p_list[1]*dh2/dy+p_list[2]*dh3/dy+p_list[3]*dh4/dy+p_list[4]*dh5/dy+p_list[5]*dh6/dy;
     193         *   dp/dz=p_list[0]*dh1/dz+p_list[1]*dh2/dz+p_list[2]*dh3/dz+p_list[3]*dh4/dz+p_list[4]*dh5/dz+p_list[5]*dh6/dz;
     194         *
     195         *   p is a vector of size 3x1 already allocated.
     196         */
     197
     198        const int NDOF3=3;
     199        const int numgrids=6;
     200        double dh1dh6[NDOF3][numgrids];
     201
     202        /*Get nodal funnctions derivatives in actual coordinate system: */
     203        GetNodalFunctionsP1Derivatives(&dh1dh6[0][0],xyz_list, gauss);
     204
     205        p[0]=this->values[0]*dh1dh6[0][0]+this->values[1]*dh1dh6[0][1]+this->values[2]*dh1dh6[0][2]+this->values[3]*dh1dh6[0][3]+this->values[4]*dh1dh6[0][4]+this->values[5]*dh1dh6[0][5];
     206        p[1]=this->values[0]*dh1dh6[1][0]+this->values[1]*dh1dh6[1][1]+this->values[2]*dh1dh6[1][2]+this->values[3]*dh1dh6[1][3]+this->values[4]*dh1dh6[1][4]+this->values[5]*dh1dh6[1][5];
     207        p[2]=this->values[0]*dh1dh6[2][0]+this->values[1]*dh1dh6[2][1]+this->values[2]*dh1dh6[2][2]+this->values[3]*dh1dh6[2][3]+this->values[4]*dh1dh6[2][4]+this->values[5]*dh1dh6[2][5];
     208
     209}
     210/*}}}*/
     211/*FUNCTION PentaVertexInput::GetVxStrainRate3d(double* epsilonvx,double* xyz_list, double* gauss,int formulation_enum) {{{1*/
     212void PentaVertexInput::GetVxStrainRate3d(double* epsilonvx,double* xyz_list, double* gauss,int formulation_enum){
     213
     214        if(formulation_enum==StokesAnalysisEnum){
     215                GetVxStrainRate3dStokes(epsilonvx,xyz_list,gauss);
     216        }
     217        else if(formulation_enum==PattynFormulationEnum){
     218                GetVxStrainRate3dPattyn(epsilonvx,xyz_list,gauss);
     219        }
     220        else{
     221                ISSMERROR("Formulation enum %i (%s) not supported yet",formulation_enum,EnumAsString(formulation_enum));
     222        }
     223}
     224/*}}}*/
     225/*FUNCTION PentaVertexInput::GetVxStrainRate3dStokes(double* epsilonvx,double* xyz_list, double* gauss) {{{1*/
     226void PentaVertexInput::GetVxStrainRate3dStokes(double* epsilonvx,double* xyz_list, double* gauss){
     227        int i,j;
     228
     229        const int numgrids=6;
     230        const int DOFVELOCITY=3;
     231        double B[8][27];
     232        double B_reduced[6][DOFVELOCITY*numgrids];
     233        double velocity[3][DOFVELOCITY];
     234
     235        /*Get B matrix: */
     236        GetBStokes(&B[0][0], xyz_list, gauss);
     237        /*Create a reduced matrix of B to get rid of pressure */
     238        for (i=0;i<6;i++){
     239                for (j=0;j<3;j++){
     240                        B_reduced[i][j]=B[i][j];
     241                }
     242                for (j=4;j<7;j++){
     243                        B_reduced[i][j-1]=B[i][j];
     244                }
     245                for (j=8;j<11;j++){
     246                        B_reduced[i][j-2]=B[i][j];
     247                }
     248                for (j=12;j<15;j++){
     249                        B_reduced[i][j-3]=B[i][j];
     250                }
     251                for (j=16;j<19;j++){
     252                        B_reduced[i][j-4]=B[i][j];
     253                }
     254                for (j=20;j<23;j++){
     255                        B_reduced[i][j-5]=B[i][j];
     256                }
     257        }
     258
     259        /*Here, we are computing the strain rate of (vx,0,0)*/
     260        for(i=0;i<numgrids;i++){
     261                velocity[i][0]=this->values[i];
     262                velocity[i][1]=0.0;
     263                velocity[i][2]=0.0;
     264        }
     265        /*Multiply B by velocity, to get strain rate: */
     266        MatrixMultiply(&B_reduced[0][0],6,DOFVELOCITY*numgrids,0,&velocity[0][0],DOFVELOCITY*numgrids,1,0,epsilonvx,0);
     267
     268}
     269/*}}}*/
     270/*FUNCTION PentaVertexInput::GetVxStrainRate3dPattyn(double* epsilonvx,double* xyz_list, double* gauss) {{{1*/
     271void PentaVertexInput::GetVxStrainRate3dPattyn(double* epsilonvx,double* xyz_list, double* gauss){
     272
     273        int i;
     274        const int numgrids=6;
     275        const int NDOF2=2;
     276        double B[5][NDOF2*numgrids];
     277        double velocity[numgrids][NDOF2];
     278
     279        /*Get B matrix: */
     280        GetBPattyn(&B[0][0], xyz_list, gauss);
     281
     282        /*Here, we are computing the strain rate of (vx,0)*/
     283        for(i=0;i<numgrids;i++){
     284                velocity[i][0]=this->values[i];
     285                velocity[i][1]=0.0;
     286        }
     287
     288        /*Multiply B by velocity, to get strain rate: */
     289        MatrixMultiply( &B[0][0],5,NDOF2*numgrids,0,
     290                                &velocity[0][0],NDOF2*numgrids,1,0,
     291                                epsilonvx,0);
     292
     293}
     294/*}}}*/
     295/*FUNCTION PentaVertexInput::GetVyStrainRate3d(double* epsilonvy,double* xyz_list, double* gauss,int formulation_enum) {{{1*/
     296void PentaVertexInput::GetVyStrainRate3d(double* epsilonvy,double* xyz_list, double* gauss,int formulation_enum){
     297
     298        if(formulation_enum==StokesAnalysisEnum){
     299                GetVyStrainRate3dStokes(epsilonvy,xyz_list,gauss);
     300        }
     301        else if(formulation_enum==PattynFormulationEnum){
     302                GetVyStrainRate3dPattyn(epsilonvy,xyz_list,gauss);
     303        }
     304        else{
     305                ISSMERROR("Formulation enum %i (%s) not supported yet",formulation_enum,EnumAsString(formulation_enum));
     306        }
     307}
     308/*}}}*/
     309/*FUNCTION PentaVertexInput::GetVyStrainRate3dStokes(double* epsilonvy,double* xyz_list, double* gauss) {{{1*/
     310void PentaVertexInput::GetVyStrainRate3dStokes(double* epsilonvy,double* xyz_list, double* gauss){
     311        int i,j;
     312
     313        const int numgrids=6;
     314        const int DOFVELOCITY=3;
     315        double B[8][27];
     316        double B_reduced[6][DOFVELOCITY*numgrids];
     317        double velocity[3][DOFVELOCITY];
     318
     319        /*Get B matrix: */
     320        GetBStokes(&B[0][0], xyz_list, gauss);
     321        /*Create a reduced matrix of B to get rid of pressure */
     322        for (i=0;i<6;i++){
     323                for (j=0;j<3;j++){
     324                        B_reduced[i][j]=B[i][j];
     325                }
     326                for (j=4;j<7;j++){
     327                        B_reduced[i][j-1]=B[i][j];
     328                }
     329                for (j=8;j<11;j++){
     330                        B_reduced[i][j-2]=B[i][j];
     331                }
     332                for (j=12;j<15;j++){
     333                        B_reduced[i][j-3]=B[i][j];
     334                }
     335                for (j=16;j<19;j++){
     336                        B_reduced[i][j-4]=B[i][j];
     337                }
     338                for (j=20;j<23;j++){
     339                        B_reduced[i][j-5]=B[i][j];
     340                }
     341        }
     342
     343        /*Here, we are computing the strain rate of (0,vy,0)*/
     344        for(i=0;i<numgrids;i++){
     345                velocity[i][0]=0.0;
     346                velocity[i][1]=this->values[i];
     347                velocity[i][2]=0.0;
     348        }
     349        /*Multiply B by velocity, to get strain rate: */
     350        MatrixMultiply(&B_reduced[0][0],6,DOFVELOCITY*numgrids,0,&velocity[0][0],DOFVELOCITY*numgrids,1,0,epsilonvy,0);
     351
     352}
     353/*}}}*/
     354/*FUNCTION PentaVertexInput::GetVyStrainRate3dPattyn(double* epsilonvy,double* xyz_list, double* gauss) {{{1*/
     355void PentaVertexInput::GetVyStrainRate3dPattyn(double* epsilonvy,double* xyz_list, double* gauss){
     356
     357        int i;
     358        const int numgrids=6;
     359        const int NDOF2=2;
     360        double B[5][NDOF2*numgrids];
     361        double velocity[numgrids][NDOF2];
     362
     363        /*Get B matrix: */
     364        GetBPattyn(&B[0][0], xyz_list, gauss);
     365
     366        /*Here, we are computing the strain rate of (0,vy)*/
     367        for(i=0;i<numgrids;i++){
     368                velocity[i][0]=0.0;
     369                velocity[i][1]=this->values[i];
     370        }
     371
     372        /*Multiply B by velocity, to get strain rate: */
     373        MatrixMultiply( &B[0][0],5,NDOF2*numgrids,0,
     374                                &velocity[0][0],NDOF2*numgrids,1,0,
     375                                epsilonvy,0);
     376
     377}
     378/*}}}*/
     379/*FUNCTION PentaVertexInput::GetVzStrainRate3d(double* epsilonvz,double* xyz_list, double* gauss,int formulation_enum) {{{1*/
     380void PentaVertexInput::GetVzStrainRate3d(double* epsilonvz,double* xyz_list, double* gauss,int formulation_enum){
     381
     382        if(formulation_enum==StokesAnalysisEnum){
     383                GetVzStrainRate3dStokes(epsilonvz,xyz_list,gauss);
     384        }
     385        else if(formulation_enum==PattynFormulationEnum){
     386                GetVzStrainRate3dPattyn(epsilonvz,xyz_list,gauss);
     387        }
     388        else{
     389                ISSMERROR("Formulation enum %i (%s) not supported yet",formulation_enum,EnumAsString(formulation_enum));
     390        }
     391}
     392/*}}}*/
     393/*FUNCTION PentaVertexInput::GetVzStrainRate3dStokes(double* epsilonvz,double* xyz_list, double* gauss) {{{1*/
     394void PentaVertexInput::GetVzStrainRate3dStokes(double* epsilonvz,double* xyz_list, double* gauss){
     395        int i,j;
     396
     397        const int numgrids=6;
     398        const int DOFVELOCITY=3;
     399        double B[8][27];
     400        double B_reduced[6][DOFVELOCITY*numgrids];
     401        double velocity[3][DOFVELOCITY];
     402
     403        /*Get B matrix: */
     404        GetBStokes(&B[0][0], xyz_list, gauss);
     405        /*Create a reduced matrix of B to get rid of pressure */
     406        for (i=0;i<6;i++){
     407                for (j=0;j<3;j++){
     408                        B_reduced[i][j]=B[i][j];
     409                }
     410                for (j=4;j<7;j++){
     411                        B_reduced[i][j-1]=B[i][j];
     412                }
     413                for (j=8;j<11;j++){
     414                        B_reduced[i][j-2]=B[i][j];
     415                }
     416                for (j=12;j<15;j++){
     417                        B_reduced[i][j-3]=B[i][j];
     418                }
     419                for (j=16;j<19;j++){
     420                        B_reduced[i][j-4]=B[i][j];
     421                }
     422                for (j=20;j<23;j++){
     423                        B_reduced[i][j-5]=B[i][j];
     424                }
     425        }
     426
     427        /*Here, we are computing the strain rate of (0,0,vz)*/
     428        for(i=0;i<numgrids;i++){
     429                velocity[i][0]=0.0;
     430                velocity[i][1]=0.0;
     431                velocity[i][2]=this->values[i];
     432        }
     433
     434        /*Multiply B by velocity, to get strain rate: */
     435        MatrixMultiply(&B_reduced[0][0],6,DOFVELOCITY*numgrids,0,&velocity[0][0],DOFVELOCITY*numgrids,1,0,epsilonvz,0);
     436
     437}
     438/*}}}*/
     439/*FUNCTION PentaVertexInput::GetVzStrainRate3dPattyn(double* epsilonvz,double* xyz_list, double* gauss) {{{1*/
     440void PentaVertexInput::GetVzStrainRate3dPattyn(double* epsilonvz,double* xyz_list, double* gauss){
     441
     442        /*vz does not contribute to the strain rate in Pattyn/Blatter's model*/
     443        for (int i=0;i<5;i++){
     444                epsilonvz[i]=0.0;
     445        }
     446
     447}
    164448/*}}}*/
    165449/*FUNCTION PentaVertexInput::ChangeEnum(int newenumtype){{{1*/
    … …  
    173457}
    174458/*}}}*/
     459
     460/*Intermediary*/
     461/*FUNCTION PentaVertexInput::GetNodalFunctionsP1 {{{1*/
     462void PentaVertexInput::GetNodalFunctionsP1(double* l1l6, double* gauss_coord){
     463
     464        /*This routine returns the values of the nodal functions  at the gaussian point.*/
     465
     466        l1l6[0]=gauss_coord[0]*(1-gauss_coord[3])/2.0;
     467
     468        l1l6[1]=gauss_coord[1]*(1-gauss_coord[3])/2.0;
     469
     470        l1l6[2]=gauss_coord[2]*(1-gauss_coord[3])/2.0;
     471
     472        l1l6[3]=gauss_coord[0]*(1+gauss_coord[3])/2.0;
     473
     474        l1l6[4]=gauss_coord[1]*(1+gauss_coord[3])/2.0;
     475
     476        l1l6[5]=gauss_coord[2]*(1+gauss_coord[3])/2.0;
     477
     478}
     479/*}}}*/
     480/*FUNCTION PentaVertexInput::GetNodalFunctionsMINI{{{1*/
     481void PentaVertexInput::GetNodalFunctionsMINI(double* l1l7, double* gauss_coord){
     482
     483        /*This routine returns the values of the nodal functions  at the gaussian point.*/
     484
     485        /*First nodal function: */
     486        l1l7[0]=gauss_coord[0]*(1.0-gauss_coord[3])/2.0;
     487
     488        /*Second nodal function: */
     489        l1l7[1]=gauss_coord[1]*(1.0-gauss_coord[3])/2.0;
     490
     491        /*Third nodal function: */
     492        l1l7[2]=gauss_coord[2]*(1.0-gauss_coord[3])/2.0;
     493
     494        /*Fourth nodal function: */
     495        l1l7[3]=gauss_coord[0]*(1.0+gauss_coord[3])/2.0;
     496
     497        /*Fifth nodal function: */
     498        l1l7[4]=gauss_coord[1]*(1.0+gauss_coord[3])/2.0;
     499
     500        /*Sixth nodal function: */
     501        l1l7[5]=gauss_coord[2]*(1.0+gauss_coord[3])/2.0;
     502
     503        /*Seventh nodal function: */
     504        l1l7[6]=27*gauss_coord[0]*gauss_coord[1]*gauss_coord[2]*(1.0+gauss_coord[3])*(1.0-gauss_coord[3]);
     505
     506}
     507/*}}}*/
     508/*FUNCTION PentaVertexInput::GetNodalFunctionsP1Derivatives {{{1*/
     509void PentaVertexInput::GetNodalFunctionsP1Derivatives(double* dh1dh6,double* xyz_list, double* gauss_coord){
     510
     511        /*This routine returns the values of the nodal functions derivatives  (with respect to the actual coordinate system: */
     512        int i;
     513        const int NDOF3=3;
     514        const int numgrids=6;
     515
     516        double dh1dh6_ref[NDOF3][numgrids];
     517        double Jinv[NDOF3][NDOF3];
     518
     519        /*Get derivative values with respect to parametric coordinate system: */
     520        GetNodalFunctionsP1DerivativesReference(&dh1dh6_ref[0][0], gauss_coord);
     521
     522        /*Get Jacobian invert: */
     523        GetJacobianInvert(&Jinv[0][0], xyz_list, gauss_coord);
     524
     525        /*Build dh1dh3:
     526         *
     527         * [dhi/dx]= Jinv*[dhi/dr]
     528         * [dhi/dy]       [dhi/ds]
     529         * [dhi/dz]       [dhi/dn]
     530         */
     531
     532        for (i=0;i<numgrids;i++){
     533                *(dh1dh6+numgrids*0+i)=Jinv[0][0]*dh1dh6_ref[0][i]+Jinv[0][1]*dh1dh6_ref[1][i]+Jinv[0][2]*dh1dh6_ref[2][i];
     534                *(dh1dh6+numgrids*1+i)=Jinv[1][0]*dh1dh6_ref[0][i]+Jinv[1][1]*dh1dh6_ref[1][i]+Jinv[1][2]*dh1dh6_ref[2][i];
     535                *(dh1dh6+numgrids*2+i)=Jinv[2][0]*dh1dh6_ref[0][i]+Jinv[2][1]*dh1dh6_ref[1][i]+Jinv[2][2]*dh1dh6_ref[2][i];
     536        }
     537
     538}
     539/*}}}*/
     540/*FUNCTION PentaVertexInput::GetNodalFunctionsMINIDerivatives{{{1*/
     541void PentaVertexInput::GetNodalFunctionsMINIDerivatives(double* dh1dh7,double* xyz_list, double* gauss_coord){
     542
     543        /*This routine returns the values of the nodal functions derivatives  (with respect to the
     544         * actual coordinate system: */
     545
     546        int i;
     547
     548        const  int numgrids=7;
     549        double dh1dh7_ref[3][numgrids];
     550        double Jinv[3][3];
     551
     552
     553        /*Get derivative values with respect to parametric coordinate system: */
     554        GetNodalFunctionsMINIDerivativesReference(&dh1dh7_ref[0][0], gauss_coord);
     555
     556        /*Get Jacobian invert: */
     557        GetJacobianInvert(&Jinv[0][0], xyz_list, gauss_coord);
     558
     559        /*Build dh1dh6:
     560         *
     561         * [dhi/dx]= Jinv'*[dhi/dr]
     562         * [dhi/dy]        [dhi/ds]
     563         * [dhi/dz]        [dhi/dzeta]
     564         */
     565
     566        for (i=0;i<numgrids;i++){
     567                *(dh1dh7+numgrids*0+i)=Jinv[0][0]*dh1dh7_ref[0][i]+Jinv[0][1]*dh1dh7_ref[1][i]+Jinv[0][2]*dh1dh7_ref[2][i];
     568                *(dh1dh7+numgrids*1+i)=Jinv[1][0]*dh1dh7_ref[0][i]+Jinv[1][1]*dh1dh7_ref[1][i]+Jinv[1][2]*dh1dh7_ref[2][i];
     569                *(dh1dh7+numgrids*2+i)=Jinv[2][0]*dh1dh7_ref[0][i]+Jinv[2][1]*dh1dh7_ref[1][i]+Jinv[2][2]*dh1dh7_ref[2][i];
     570        }
     571
     572}
     573/*}}}*/
     574/*FUNCTION PentaVertexInput::GetNodalFunctionsP1DerivativesReference {{{1*/
     575void PentaVertexInput::GetNodalFunctionsP1DerivativesReference(double* dl1dl6,double* gauss_coord){
     576
     577        /*This routine returns the values of the nodal functions derivatives  (with respect to the
     578         * natural coordinate system) at the gaussian point. Those values vary along xi,eta,z */
     579
     580        const int numgrids=6;
     581        double A1,A2,A3,z;
     582
     583        A1=gauss_coord[0]; //first area coordinate value. In term of xi and eta: A1=(1-xi)/2-eta/(2*SQRT3);
     584        A2=gauss_coord[1]; //second area coordinate value In term of xi and eta: A2=(1+xi)/2-eta/(2*SQRT3);
     585        A3=gauss_coord[2]; //third area coordinate value  In term of xi and eta: A3=y/SQRT3;
     586        z=gauss_coord[3]; //fourth vertical coordinate value. Corresponding nodal function: (1-z)/2 and (1+z)/2
     587
     588
     589        /*First nodal function derivatives. The corresponding nodal function is N=A1*(1-z)/2. Its derivatives follow*/
     590        *(dl1dl6+numgrids*0+0)=-0.5*(1.0-z)/2.0;
     591        *(dl1dl6+numgrids*1+0)=-0.5/SQRT3*(1.0-z)/2.0;
     592        *(dl1dl6+numgrids*2+0)=-0.5*A1;
     593
     594        /*Second nodal function: The corresponding nodal function is N=A2*(1-z)/2. Its derivatives follow*/
     595        *(dl1dl6+numgrids*0+1)=0.5*(1.0-z)/2.0;
     596        *(dl1dl6+numgrids*1+1)=-0.5/SQRT3*(1.0-z)/2.0;
     597        *(dl1dl6+numgrids*2+1)=-0.5*A2;
     598
     599        /*Third nodal function: The corresponding nodal function is N=A3*(1-z)/2. Its derivatives follow*/
     600        *(dl1dl6+numgrids*0+2)=0.0;
     601        *(dl1dl6+numgrids*1+2)=1.0/SQRT3*(1.0-z)/2.0;
     602        *(dl1dl6+numgrids*2+2)=-0.5*A3;
     603
     604        /*Fourth nodal function: The corresponding nodal function is N=A1*(1+z)/2. Its derivatives follow*/
     605        *(dl1dl6+numgrids*0+3)=-0.5*(1.0+z)/2.0;
     606        *(dl1dl6+numgrids*1+3)=-0.5/SQRT3*(1.0+z)/2.0;
     607        *(dl1dl6+numgrids*2+3)=0.5*A1;
     608
     609        /*Fifth nodal function: The corresponding nodal function is N=A2*(1+z)/2. Its derivatives follow*/
     610        *(dl1dl6+numgrids*0+4)=0.5*(1.0+z)/2.0;
     611        *(dl1dl6+numgrids*1+4)=-0.5/SQRT3*(1.0+z)/2.0;
     612        *(dl1dl6+numgrids*2+4)=0.5*A2;
     613
     614        /*Sixth nodal function: The corresponding nodal function is N=A3*(1+z)/2. Its derivatives follow*/
     615        *(dl1dl6+numgrids*0+5)=0.0;
     616        *(dl1dl6+numgrids*1+5)=1.0/SQRT3*(1.0+z)/2.0;
     617        *(dl1dl6+numgrids*2+5)=0.5*A3;
     618}
     619/*}}}*/
     620/*FUNCTION PentaVertexInput::GetNodalFunctionsMINIDerivativesReference{{{1*/
     621void PentaVertexInput::GetNodalFunctionsMINIDerivativesReference(double* dl1dl7,double* gauss_coord){
     622
     623        /*This routine returns the values of the nodal functions derivatives  (with respect to the
     624         * natural coordinate system) at the gaussian point. */
     625
     626        int    numgrids=7; //six plus bubble grids
     627
     628        double r=gauss_coord[1]-gauss_coord[0];
     629        double s=-3.0/SQRT3*(gauss_coord[0]+gauss_coord[1]-2.0/3.0);
     630        double zeta=gauss_coord[3];
     631
     632        /*First nodal function: */
     633        *(dl1dl7+numgrids*0+0)=-0.5*(1.0-zeta)/2.0;
     634        *(dl1dl7+numgrids*1+0)=-SQRT3/6.0*(1.0-zeta)/2.0;
     635        *(dl1dl7+numgrids*2+0)=-0.5*(-0.5*r-SQRT3/6.0*s+ONETHIRD);
     636
     637        /*Second nodal function: */
     638        *(dl1dl7+numgrids*0+1)=0.5*(1.0-zeta)/2.0;
     639        *(dl1dl7+numgrids*1+1)=-SQRT3/6.0*(1.0-zeta)/2.0;
     640        *(dl1dl7+numgrids*2+1)=-0.5*(0.5*r-SQRT3/6.0*s+ONETHIRD);
     641
     642        /*Third nodal function: */
     643        *(dl1dl7+numgrids*0+2)=0;
     644        *(dl1dl7+numgrids*1+2)=SQRT3/3.0*(1.0-zeta)/2.0;
     645        *(dl1dl7+numgrids*2+2)=-0.5*(SQRT3/3.0*s+ONETHIRD);
     646
     647        /*Fourth nodal function: */
     648        *(dl1dl7+numgrids*0+3)=-0.5*(1.0+zeta)/2.0;
     649        *(dl1dl7+numgrids*1+3)=-SQRT3/6.0*(1.0+zeta)/2.0;
     650        *(dl1dl7+numgrids*2+3)=0.5*(-0.5*r-SQRT3/6.0*s+ONETHIRD);
     651
     652        /*Fith nodal function: */
     653        *(dl1dl7+numgrids*0+4)=0.5*(1.0+zeta)/2.0;
     654        *(dl1dl7+numgrids*1+4)=-SQRT3/6.0*(1.0+zeta)/2.0;
     655        *(dl1dl7+numgrids*2+4)=0.5*(0.5*r-SQRT3/6.0*s+ONETHIRD);
     656
     657        /*Sixth nodal function: */
     658        *(dl1dl7+numgrids*0+5)=0;
     659        *(dl1dl7+numgrids*1+5)=SQRT3/3.0*(1.0+zeta)/2.0;
     660        *(dl1dl7+numgrids*2+5)=0.5*(SQRT3/3.0*s+ONETHIRD);
     661
     662        /*Seventh nodal function: */
     663        *(dl1dl7+numgrids*0+6)=9.0/2.0*r*(1.0+zeta)*(zeta-1.0)*(SQRT3*s+1.0);
     664        *(dl1dl7+numgrids*1+6)=9.0/4.0*(1+zeta)*(1-zeta)*(SQRT3*pow(s,2.0)-2.0*s-SQRT3*pow(r,2.0));
     665        *(dl1dl7+numgrids*2+6)=27*gauss_coord[0]*gauss_coord[1]*gauss_coord[2]*(-2.0*zeta);
     666
     667}
     668/*}}}*/
     669/*FUNCTION PentaVertexInput::GetJacobian {{{1*/
     670void PentaVertexInput::GetJacobian(double* J, double* xyz_list,double* gauss_coord){
     671
     672        const int NDOF3=3;
     673        int i,j;
     674
     675        /*The Jacobian is constant over the element, discard the gaussian points.
     676         * J is assumed to have been allocated of size NDOF2xNDOF2.*/
     677
     678        double A1,A2,A3; //area coordinates
     679        double xi,eta,zi; //parametric coordinates
     680
     681        double x1,x2,x3,x4,x5,x6;
     682        double y1,y2,y3,y4,y5,y6;
     683        double z1,z2,z3,z4,z5,z6;
     684
     685        /*Figure out xi,eta and zi (parametric coordinates), for this gaussian point: */
     686        A1=gauss_coord[0];
     687        A2=gauss_coord[1];
     688        A3=gauss_coord[2];
     689
     690        xi=A2-A1;
     691        eta=SQRT3*A3;
     692        zi=gauss_coord[3];
     693
     694        x1=*(xyz_list+3*0+0);
     695        x2=*(xyz_list+3*1+0);
     696        x3=*(xyz_list+3*2+0);
     697        x4=*(xyz_list+3*3+0);
     698        x5=*(xyz_list+3*4+0);
     699        x6=*(xyz_list+3*5+0);
     700
     701        y1=*(xyz_list+3*0+1);
     702        y2=*(xyz_list+3*1+1);
     703        y3=*(xyz_list+3*2+1);
     704        y4=*(xyz_list+3*3+1);
     705        y5=*(xyz_list+3*4+1);
     706        y6=*(xyz_list+3*5+1);
     707
     708        z1=*(xyz_list+3*0+2);
     709        z2=*(xyz_list+3*1+2);
     710        z3=*(xyz_list+3*2+2);
     711        z4=*(xyz_list+3*3+2);
     712        z5=*(xyz_list+3*4+2);
     713        z6=*(xyz_list+3*5+2);
     714
     715        *(J+NDOF3*0+0)=0.25*(x1-x2-x4+x5)*zi+0.25*(-x1+x2-x4+x5);
     716        *(J+NDOF3*1+0)=SQRT3/12.0*(x1+x2-2*x3-x4-x5+2*x6)*zi+SQRT3/12.0*(-x1-x2+2*x3-x4-x5+2*x6);
     717        *(J+NDOF3*2+0)=SQRT3/12.0*(x1+x2-2*x3-x4-x5+2*x6)*eta+1/4*(x1-x2-x4+x5)*xi +0.25*(-x1+x5-x2+x4);
     718
     719        *(J+NDOF3*0+1)=0.25*(y1-y2-y4+y5)*zi+0.25*(-y1+y2-y4+y5);
     720        *(J+NDOF3*1+1)=SQRT3/12.0*(y1+y2-2*y3-y4-y5+2*y6)*zi+SQRT3/12.0*(-y1-y2+2*y3-y4-y5+2*y6);
     721        *(J+NDOF3*2+1)=SQRT3/12.0*(y1+y2-2*y3-y4-y5+2*y6)*eta+0.25*(y1-y2-y4+y5)*xi+0.25*(y4-y1+y5-y2);
     722
     723        *(J+NDOF3*0+2)=0.25*(z1-z2-z4+z5)*zi+0.25*(-z1+z2-z4+z5);
     724        *(J+NDOF3*1+2)=SQRT3/12.0*(z1+z2-2*z3-z4-z5+2*z6)*zi+SQRT3/12.0*(-z1-z2+2*z3-z4-z5+2*z6);
     725        *(J+NDOF3*2+2)=SQRT3/12.0*(z1+z2-2*z3-z4-z5+2*z6)*eta+0.25*(z1-z2-z4+z5)*xi+0.25*(-z1+z5-z2+z4);
     726
     727}
     728/*}}}*/
     729/*FUNCTION PentaVertexInput::GetJacobianInvert {{{1*/
     730void PentaVertexInput::GetJacobianInvert(double*  Jinv, double* xyz_list,double* gauss_coord){
     731
     732        double Jdet;
     733        const int NDOF3=3;
     734
     735        /*Call Jacobian routine to get the jacobian:*/
     736        GetJacobian(Jinv, xyz_list, gauss_coord);
     737
     738        /*Invert Jacobian matrix: */
     739        MatrixInverse(Jinv,NDOF3,NDOF3,NULL,0,&Jdet);
     740}
     741/*}}}*/
     742/*FUNCTION PentaVertexInput::GetBPattyn {{{1*/
     743void PentaVertexInput::GetBPattyn(double* B, double* xyz_list, double* gauss_coord){
     744        /*Compute B  matrix. B=[B1 B2 B3 B4 B5 B6] where Bi is of size 5*NDOF2.
     745         * For grid i, Bi can be expressed in the actual coordinate system
     746         * by:
     747         *       Bi=[ dh/dx          0      ]
     748         *          [   0           dh/dy   ]
     749         *          [ 1/2*dh/dy  1/2*dh/dx  ]
     750         *          [ 1/2*dh/dz      0      ]
     751         *          [  0         1/2*dh/dz  ]
     752         * where h is the interpolation function for grid i.
     753         *
     754         * We assume B has been allocated already, of size: 5x(NDOF2*numgrids)
     755         */
     756
     757        int i;
     758        const int numgrids=6;
     759        const int NDOF3=3;
     760        const int NDOF2=2;
     761
     762        double dh1dh6[NDOF3][numgrids];
     763
     764        /*Get dh1dh6 in actual coordinate system: */
     765        GetNodalFunctionsP1Derivatives(&dh1dh6[0][0],xyz_list, gauss_coord);
     766
     767        /*Build B: */
     768        for (i=0;i<numgrids;i++){
     769                *(B+NDOF2*numgrids*0+NDOF2*i)=dh1dh6[0][i];
     770                *(B+NDOF2*numgrids*0+NDOF2*i+1)=0.0;
     771
     772                *(B+NDOF2*numgrids*1+NDOF2*i)=0.0;
     773                *(B+NDOF2*numgrids*1+NDOF2*i+1)=dh1dh6[1][i];
     774
     775                *(B+NDOF2*numgrids*2+NDOF2*i)=(float).5*dh1dh6[1][i];
     776                *(B+NDOF2*numgrids*2+NDOF2*i+1)=(float).5*dh1dh6[0][i];
     777
     778                *(B+NDOF2*numgrids*3+NDOF2*i)=(float).5*dh1dh6[2][i];
     779                *(B+NDOF2*numgrids*3+NDOF2*i+1)=0.0;
     780
     781                *(B+NDOF2*numgrids*4+NDOF2*i)=0.0;
     782                *(B+NDOF2*numgrids*4+NDOF2*i+1)=(float).5*dh1dh6[2][i];
     783        }
     784
     785}
     786/*}}}*/
     787/*FUNCTION PentaVertexInput::GetBStokes {{{1*/
     788void PentaVertexInput::GetBStokes(double* B, double* xyz_list, double* gauss_coord){
     789
     790        /*Compute B  matrix. B=[B1 B2 B3 B4 B5 B6] where Bi is of size 3*DOFPERGRID.
     791         * For grid i, Bi can be expressed in the actual coordinate system
     792         * by:          Bi=[ dh/dx          0             0       0  ]
     793         *                                      [   0           dh/dy           0       0  ]
     794         *                                      [   0             0           dh/dy     0  ]
     795         *                                      [ 1/2*dh/dy    1/2*dh/dx        0       0  ]
     796         *                                      [ 1/2*dh/dz       0         1/2*dh/dx   0  ]
     797         *                                      [   0          1/2*dh/dz    1/2*dh/dy   0  ]
     798         *                                      [   0             0             0       h  ]
     799         *                                      [ dh/dx         dh/dy         dh/dz     0  ]
     800         *      where h is the interpolation function for grid i.
     801         *      Same thing for Bb except the last column that does not exist.
     802         */
     803
     804        int i;
     805        const int calculationdof=3;
     806        const int numgrids=6;
     807        int DOFPERGRID=4;
     808
     809        double dh1dh7[calculationdof][numgrids+1];
     810        double l1l6[numgrids];
     811
     812
     813        /*Get dh1dh7 in actual coordinate system: */
     814        GetNodalFunctionsMINIDerivatives(&dh1dh7[0][0],xyz_list, gauss_coord);
     815        GetNodalFunctionsP1(l1l6, gauss_coord);
     816
     817        /*Build B: */
     818        for (i=0;i<numgrids+1;i++){
     819                *(B+(DOFPERGRID*numgrids+3)*0+DOFPERGRID*i)=dh1dh7[0][i]; //B[0][DOFPERGRID*i]=dh1dh6[0][i];
     820                *(B+(DOFPERGRID*numgrids+3)*0+DOFPERGRID*i+1)=0;
     821                *(B+(DOFPERGRID*numgrids+3)*0+DOFPERGRID*i+2)=0;
     822                *(B+(DOFPERGRID*numgrids+3)*1+DOFPERGRID*i)=0;
     823                *(B+(DOFPERGRID*numgrids+3)*1+DOFPERGRID*i+1)=dh1dh7[1][i];
     824                *(B+(DOFPERGRID*numgrids+3)*1+DOFPERGRID*i+2)=0;
     825                *(B+(DOFPERGRID*numgrids+3)*2+DOFPERGRID*i)=0;
     826                *(B+(DOFPERGRID*numgrids+3)*2+DOFPERGRID*i+1)=0;
     827                *(B+(DOFPERGRID*numgrids+3)*2+DOFPERGRID*i+2)=dh1dh7[2][i];
     828                *(B+(DOFPERGRID*numgrids+3)*3+DOFPERGRID*i)=(float).5*dh1dh7[1][i];
     829                *(B+(DOFPERGRID*numgrids+3)*3+DOFPERGRID*i+1)=(float).5*dh1dh7[0][i];
     830                *(B+(DOFPERGRID*numgrids+3)*3+DOFPERGRID*i+2)=0;
     831                *(B+(DOFPERGRID*numgrids+3)*4+DOFPERGRID*i)=(float).5*dh1dh7[2][i];
     832                *(B+(DOFPERGRID*numgrids+3)*4+DOFPERGRID*i+1)=0;
     833                *(B+(DOFPERGRID*numgrids+3)*4+DOFPERGRID*i+2)=(float).5*dh1dh7[0][i];
     834                *(B+(DOFPERGRID*numgrids+3)*5+DOFPERGRID*i)=0;
     835                *(B+(DOFPERGRID*numgrids+3)*5+DOFPERGRID*i+1)=(float).5*dh1dh7[2][i];
     836                *(B+(DOFPERGRID*numgrids+3)*5+DOFPERGRID*i+2)=(float).5*dh1dh7[1][i];
     837                *(B+(DOFPERGRID*numgrids+3)*6+DOFPERGRID*i)=0;
     838                *(B+(DOFPERGRID*numgrids+3)*6+DOFPERGRID*i+1)=0;
     839                *(B+(DOFPERGRID*numgrids+3)*6+DOFPERGRID*i+2)=0;
     840                *(B+(DOFPERGRID*numgrids+3)*7+DOFPERGRID*i)=dh1dh7[0][i];
     841                *(B+(DOFPERGRID*numgrids+3)*7+DOFPERGRID*i+1)=dh1dh7[1][i];
     842                *(B+(DOFPERGRID*numgrids+3)*7+DOFPERGRID*i+2)=dh1dh7[2][i];
     843        }
     844
     845        for (i=0;i<numgrids;i++){ //last column not for the bubble function
     846                *(B+(DOFPERGRID*numgrids+3)*0+DOFPERGRID*i+3)=0;
     847                *(B+(DOFPERGRID*numgrids+3)*1+DOFPERGRID*i+3)=0;
     848                *(B+(DOFPERGRID*numgrids+3)*2+DOFPERGRID*i+3)=0;
     849                *(B+(DOFPERGRID*numgrids+3)*3+DOFPERGRID*i+3)=0;
     850                *(B+(DOFPERGRID*numgrids+3)*4+DOFPERGRID*i+3)=0;
     851                *(B+(DOFPERGRID*numgrids+3)*5+DOFPERGRID*i+3)=0;
     852                *(B+(DOFPERGRID*numgrids+3)*6+DOFPERGRID*i+3)=l1l6[i];
     853                *(B+(DOFPERGRID*numgrids+3)*7+DOFPERGRID*i+3)=0;
     854        }
     855
     856}
     857/*}}}*/

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

    r3830 r3840  
    6464                void GetVxStrainRate2d(double* epsilonvx,double* xyz_list, double* gauss,int formulation_enum=MacAyealFormulationEnum){ISSMERROR("not implemented yet");};
    6565                void GetVyStrainRate2d(double* epsilonvy,double* xyz_list, double* gauss,int formulation_enum=MacAyealFormulationEnum){ISSMERROR("not implemented yet");};
    66                 void GetVxStrainRate3d(double* epsilonvx,double* xyz_list, double* gauss,int formulation_enum=StokesFormulationEnum){ISSMERROR("not implemented yet");};
    67                 void GetVyStrainRate3d(double* epsilonvy,double* xyz_list, double* gauss,int formulation_enum=StokesFormulationEnum){ISSMERROR("not implemented yet");};
    68                 void GetVzStrainRate3d(double* epsilonvz,double* xyz_list, double* gauss,int formulation_enum=StokesFormulationEnum){ISSMERROR("not implemented yet");};
     66                void GetVxStrainRate3d(double* epsilonvx,double* xyz_list, double* gauss,int formulation_enum=StokesFormulationEnum);
     67                void GetVxStrainRate3dStokes(double* epsilonvx,double* xyz_list, double* gauss);
     68                void GetVxStrainRate3dPattyn(double* epsilonvx,double* xyz_list, double* gauss);
     69                void GetVyStrainRate3d(double* epsilonvy,double* xyz_list, double* gauss,int formulation_enum=StokesFormulationEnum);
     70                void GetVyStrainRate3dStokes(double* epsilonvy,double* xyz_list, double* gauss);
     71                void GetVyStrainRate3dPattyn(double* epsilonvy,double* xyz_list, double* gauss);
     72                void GetVzStrainRate3d(double* epsilonvz,double* xyz_list, double* gauss,int formulation_enum=StokesFormulationEnum);
     73                void GetVzStrainRate3dStokes(double* epsilonvz,double* xyz_list, double* gauss);
     74                void GetVzStrainRate3dPattyn(double* epsilonvz,double* xyz_list, double* gauss);
    6975                void ChangeEnum(int newenumtype);
     76                void GetNodalFunctionsP1(double* l1l6, double* gauss_coord);
     77                void GetNodalFunctionsMINI(double* l1l7, double* gauss_coord);
     78                void GetNodalFunctionsP1Derivatives(double* dh1dh6,double* xyz_list, double* gauss_coord);
     79                void GetNodalFunctionsMINIDerivatives(double* dh1dh7,double* xyz_list, double* gauss_coord);
     80                void GetNodalFunctionsP1DerivativesReference(double* dl1dl6,double* gauss_coord);
     81                void GetNodalFunctionsMINIDerivativesReference(double* dl1dl7,double* gauss_coord);
     82                void GetJacobian(double* J, double* xyz_list,double* gauss_coord);
     83                void GetJacobianInvert(double*  Jinv, double* xyz_list,double* gauss_coord);
     84                void GetBPattyn(double* B, double* xyz_list, double* gauss_coord);
     85                void GetBStokes(double* B, double* xyz_list, double* gauss_coord);
    7086                /*}}}*/
    7187

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

    r3830 r3840  
    150150/*FUNCTION TriaVertexInput::GetParameterValues(double* values,double* gauss_pointers, int numgauss){{{1*/
    151151void TriaVertexInput::GetParameterValues(double* values,double* gauss_pointers, int numgauss){
    152         /*It is assume that output values has been correctly allocated*/
     152        /*It is assumed that output values has been correctly allocated*/
    153153
    154154        int i,j;
Note: See TracChangeset for help on using the changeset viewer.