Changeset 15458


Ignore:
Timestamp:
07/08/13 10:18:10 (12 years ago)
Author:
Mathieu Morlighem
Message:

CHG: express nodal function derivatives as functions of area coordinates

File:
1 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/classes/Elements/PentaRef.cpp

    r15436 r15458  
    11951195
    11961196        /*Get current coordinates in reference element*/
    1197         IssmDouble r=gauss->coord2-gauss->coord1;
    1198         IssmDouble s=-3.0/SQRT3*(gauss->coord1+gauss->coord2-2.0/3.0);
    11991197        IssmDouble zeta=gauss->coord4;
    12001198
    … …  
    12041202                        dbasis[NUMNODESP1*0+0]=-0.5*(1.0-zeta)/2.0;
    12051203                        dbasis[NUMNODESP1*1+0]=-SQRT3/6.0*(1.0-zeta)/2.0;
    1206                         dbasis[NUMNODESP1*2+0]=-0.5*(-0.5*r-SQRT3/6.0*s+ONETHIRD);
     1204                        dbasis[NUMNODESP1*2+0]=-0.5*gauss->coord1;
    12071205                        /*Nodal function 2*/
    12081206                        dbasis[NUMNODESP1*0+1]=0.5*(1.0-zeta)/2.0;
    12091207                        dbasis[NUMNODESP1*1+1]=-SQRT3/6.0*(1.0-zeta)/2.0;
    1210                         dbasis[NUMNODESP1*2+1]=-0.5*(0.5*r-SQRT3/6.0*s+ONETHIRD);
     1208                        dbasis[NUMNODESP1*2+1]=-0.5*gauss->coord2;
    12111209                        /*Nodal function 3*/
    1212                         dbasis[NUMNODESP1*0+2]=0;
     1210                        dbasis[NUMNODESP1*0+2]=0.;
    12131211                        dbasis[NUMNODESP1*1+2]=SQRT3/3.0*(1.0-zeta)/2.0;
    1214                         dbasis[NUMNODESP1*2+2]=-0.5*(SQRT3/3.0*s+ONETHIRD);
     1212                        dbasis[NUMNODESP1*2+2]=-0.5*gauss->coord3;
    12151213                        /*Nodal function 4*/
    12161214                        dbasis[NUMNODESP1*0+3]=-0.5*(1.0+zeta)/2.0;
    12171215                        dbasis[NUMNODESP1*1+3]=-SQRT3/6.0*(1.0+zeta)/2.0;
    1218                         dbasis[NUMNODESP1*2+3]=0.5*(-0.5*r-SQRT3/6.0*s+ONETHIRD);
     1216                        dbasis[NUMNODESP1*2+3]=0.5*gauss->coord1;
    12191217                        /*Nodal function 5*/
    12201218                        dbasis[NUMNODESP1*0+4]=0.5*(1.0+zeta)/2.0;
    12211219                        dbasis[NUMNODESP1*1+4]=-SQRT3/6.0*(1.0+zeta)/2.0;
    1222                         dbasis[NUMNODESP1*2+4]=0.5*(0.5*r-SQRT3/6.0*s+ONETHIRD);
     1220                        dbasis[NUMNODESP1*2+4]=0.5*gauss->coord2;
    12231221                        /*Nodal function 6*/
    1224                         dbasis[NUMNODESP1*0+5]=0;
     1222                        dbasis[NUMNODESP1*0+5]=0.;
    12251223                        dbasis[NUMNODESP1*1+5]=SQRT3/3.0*(1.0+zeta)/2.0;
    1226                         dbasis[NUMNODESP1*2+5]=0.5*(SQRT3/3.0*s+ONETHIRD);
     1224                        dbasis[NUMNODESP1*2+5]=0.5*gauss->coord3;
    12271225                case MINIEnum:
    12281226                        /*Nodal function 1*/
    12291227                        dbasis[NUMNODESMINI*0+0]=-0.5*(1.0-zeta)/2.0;
    12301228                        dbasis[NUMNODESMINI*1+0]=-SQRT3/6.0*(1.0-zeta)/2.0;
    1231                         dbasis[NUMNODESMINI*2+0]=-0.5*(-0.5*r-SQRT3/6.0*s+ONETHIRD);
     1229                        dbasis[NUMNODESMINI*2+0]=-0.5*gauss->coord1;
    12321230                        /*Nodal function 2*/
    12331231                        dbasis[NUMNODESMINI*0+1]=0.5*(1.0-zeta)/2.0;
    12341232                        dbasis[NUMNODESMINI*1+1]=-SQRT3/6.0*(1.0-zeta)/2.0;
    1235                         dbasis[NUMNODESMINI*2+1]=-0.5*(0.5*r-SQRT3/6.0*s+ONETHIRD);
     1233                        dbasis[NUMNODESMINI*2+1]=-0.5*gauss->coord2;
    12361234                        /*Nodal function 3*/
    1237                         dbasis[NUMNODESMINI*0+2]=0;
     1235                        dbasis[NUMNODESMINI*0+2]=0.;
    12381236                        dbasis[NUMNODESMINI*1+2]=SQRT3/3.0*(1.0-zeta)/2.0;
    1239                         dbasis[NUMNODESMINI*2+2]=-0.5*(SQRT3/3.0*s+ONETHIRD);
     1237                        dbasis[NUMNODESMINI*2+2]=-0.5*gauss->coord3;
    12401238                        /*Nodal function 4*/
    12411239                        dbasis[NUMNODESMINI*0+3]=-0.5*(1.0+zeta)/2.0;
    12421240                        dbasis[NUMNODESMINI*1+3]=-SQRT3/6.0*(1.0+zeta)/2.0;
    1243                         dbasis[NUMNODESMINI*2+3]=0.5*(-0.5*r-SQRT3/6.0*s+ONETHIRD);
     1241                        dbasis[NUMNODESMINI*2+3]=0.5*gauss->coord1;
    12441242                        /*Nodal function 5*/
    12451243                        dbasis[NUMNODESMINI*0+4]=0.5*(1.0+zeta)/2.0;
    12461244                        dbasis[NUMNODESMINI*1+4]=-SQRT3/6.0*(1.0+zeta)/2.0;
    1247                         dbasis[NUMNODESMINI*2+4]=0.5*(0.5*r-SQRT3/6.0*s+ONETHIRD);
     1245                        dbasis[NUMNODESMINI*2+4]=0.5*gauss->coord2;
    12481246                        /*Nodal function 6*/
    1249                         dbasis[NUMNODESMINI*0+5]=0;
     1247                        dbasis[NUMNODESMINI*0+5]=0.;
    12501248                        dbasis[NUMNODESMINI*1+5]=SQRT3/3.0*(1.0+zeta)/2.0;
    1251                         dbasis[NUMNODESMINI*2+5]=0.5*(SQRT3/3.0*s+ONETHIRD);
     1249                        dbasis[NUMNODESMINI*2+5]=0.5*gauss->coord3;
    12521250                        /*Nodal function 7*/
    1253                         dbasis[NUMNODESMINI*0+6]=9.0/2.0*r*(1.0+zeta)*(zeta-1.0)*(SQRT3*s+1.0);
    1254                         dbasis[NUMNODESMINI*1+6]=9.0/4.0*(1+zeta)*(1-zeta)*(SQRT3*pow(s,2.0)-2.0*s-SQRT3*pow(r,2.0));
     1251                        dbasis[NUMNODESMINI*0+6]=27.*(1.+zeta)*(1.-zeta)*(-.5*gauss->coord2*gauss->coord3 + .5*gauss->coord1*gauss->coord3);
     1252                        dbasis[NUMNODESMINI*1+6]=27.*(1.+zeta)*(1.-zeta)*SQRT3*(-1./6.*gauss->coord2*gauss->coord3 - 1./6.*gauss->coord1*gauss->coord3 +1./3.*gauss->coord1*gauss->coord2);
    12551253                        dbasis[NUMNODESMINI*2+6]=27*gauss->coord1*gauss->coord2*gauss->coord3*(-2.0*zeta);
    12561254                        return;
    … …  
    12821280         * actual coordinate system): */
    12831281
    1284         int       i;
    12851282        IssmDouble    dh1dh7_ref[3][NUMNODESMINI];
    12861283        IssmDouble    Jinv[3][3];
    … …  
    12991296         */
    13001297
    1301         for (i=0;i<NUMNODESMINI;i++){
     1298        for(int i=0;i<NUMNODESMINI;i++){
    13021299                *(dh1dh7+NUMNODESMINI*0+i)=Jinv[0][0]*dh1dh7_ref[0][i]+Jinv[0][1]*dh1dh7_ref[1][i]+Jinv[0][2]*dh1dh7_ref[2][i];
    13031300                *(dh1dh7+NUMNODESMINI*1+i)=Jinv[1][0]*dh1dh7_ref[0][i]+Jinv[1][1]*dh1dh7_ref[1][i]+Jinv[1][2]*dh1dh7_ref[2][i];
    … …  
    13081305/*}}}*/
    13091306/*FUNCTION PentaRef::GetNodalFunctionsMINIDerivativesReference{{{*/
    1310 void PentaRef::GetNodalFunctionsMINIDerivativesReference(IssmDouble* dl1dl7,GaussPenta* gauss){
     1307void PentaRef::GetNodalFunctionsMINIDerivativesReference(IssmDouble* dbasis,GaussPenta* gauss){
    13111308
    13121309        /*This routine returns the values of the nodal functions derivatives  (with respect to the
    13131310         * natural coordinate system) at the gaussian point. */
    1314         IssmDouble r=gauss->coord2-gauss->coord1;
    1315         IssmDouble s=-3.0/SQRT3*(gauss->coord1+gauss->coord2-2.0/3.0);
    13161311        IssmDouble zeta=gauss->coord4;
    13171312
    1318         /*First nodal function: */
    1319         *(dl1dl7+NUMNODESMINI*0+0)=-0.5*(1.0-zeta)/2.0;
    1320         *(dl1dl7+NUMNODESMINI*1+0)=-SQRT3/6.0*(1.0-zeta)/2.0;
    1321         *(dl1dl7+NUMNODESMINI*2+0)=-0.5*(-0.5*r-SQRT3/6.0*s+ONETHIRD);
    1322 
    1323         /*Second nodal function: */
    1324         *(dl1dl7+NUMNODESMINI*0+1)=0.5*(1.0-zeta)/2.0;
    1325         *(dl1dl7+NUMNODESMINI*1+1)=-SQRT3/6.0*(1.0-zeta)/2.0;
    1326         *(dl1dl7+NUMNODESMINI*2+1)=-0.5*(0.5*r-SQRT3/6.0*s+ONETHIRD);
    1327 
    1328         /*Third nodal function: */
    1329         *(dl1dl7+NUMNODESMINI*0+2)=0;
    1330         *(dl1dl7+NUMNODESMINI*1+2)=SQRT3/3.0*(1.0-zeta)/2.0;
    1331         *(dl1dl7+NUMNODESMINI*2+2)=-0.5*(SQRT3/3.0*s+ONETHIRD);
    1332 
    1333         /*Fourth nodal function: */
    1334         *(dl1dl7+NUMNODESMINI*0+3)=-0.5*(1.0+zeta)/2.0;
    1335         *(dl1dl7+NUMNODESMINI*1+3)=-SQRT3/6.0*(1.0+zeta)/2.0;
    1336         *(dl1dl7+NUMNODESMINI*2+3)=0.5*(-0.5*r-SQRT3/6.0*s+ONETHIRD);
    1337 
    1338         /*Fith nodal function: */
    1339         *(dl1dl7+NUMNODESMINI*0+4)=0.5*(1.0+zeta)/2.0;
    1340         *(dl1dl7+NUMNODESMINI*1+4)=-SQRT3/6.0*(1.0+zeta)/2.0;
    1341         *(dl1dl7+NUMNODESMINI*2+4)=0.5*(0.5*r-SQRT3/6.0*s+ONETHIRD);
    1342 
    1343         /*Sixth nodal function: */
    1344         *(dl1dl7+NUMNODESMINI*0+5)=0;
    1345         *(dl1dl7+NUMNODESMINI*1+5)=SQRT3/3.0*(1.0+zeta)/2.0;
    1346         *(dl1dl7+NUMNODESMINI*2+5)=0.5*(SQRT3/3.0*s+ONETHIRD);
    1347 
    1348         /*Seventh nodal function: */
    1349         *(dl1dl7+NUMNODESMINI*0+6)=9.0/2.0*r*(1.0+zeta)*(zeta-1.0)*(SQRT3*s+1.0);
    1350         *(dl1dl7+NUMNODESMINI*1+6)=9.0/4.0*(1+zeta)*(1-zeta)*(SQRT3*pow(s,2.0)-2.0*s-SQRT3*pow(r,2.0));
    1351         *(dl1dl7+NUMNODESMINI*2+6)=27*gauss->coord1*gauss->coord2*gauss->coord3*(-2.0*zeta);
    1352 
     1313        /*Nodal function 1*/
     1314        dbasis[NUMNODESMINI*0+0]=-0.5*(1.0-zeta)/2.0;
     1315        dbasis[NUMNODESMINI*1+0]=-SQRT3/6.0*(1.0-zeta)/2.0;
     1316        dbasis[NUMNODESMINI*2+0]=-0.5*gauss->coord1;
     1317        /*Nodal function 2*/
     1318        dbasis[NUMNODESMINI*0+1]=0.5*(1.0-zeta)/2.0;
     1319        dbasis[NUMNODESMINI*1+1]=-SQRT3/6.0*(1.0-zeta)/2.0;
     1320        dbasis[NUMNODESMINI*2+1]=-0.5*gauss->coord2;
     1321        /*Nodal function 3*/
     1322        dbasis[NUMNODESMINI*0+2]=0.;
     1323        dbasis[NUMNODESMINI*1+2]=SQRT3/3.0*(1.0-zeta)/2.0;
     1324        dbasis[NUMNODESMINI*2+2]=-0.5*gauss->coord3;
     1325        /*Nodal function 4*/
     1326        dbasis[NUMNODESMINI*0+3]=-0.5*(1.0+zeta)/2.0;
     1327        dbasis[NUMNODESMINI*1+3]=-SQRT3/6.0*(1.0+zeta)/2.0;
     1328        dbasis[NUMNODESMINI*2+3]=0.5*gauss->coord1;
     1329        /*Nodal function 5*/
     1330        dbasis[NUMNODESMINI*0+4]=0.5*(1.0+zeta)/2.0;
     1331        dbasis[NUMNODESMINI*1+4]=-SQRT3/6.0*(1.0+zeta)/2.0;
     1332        dbasis[NUMNODESMINI*2+4]=0.5*gauss->coord2;
     1333        /*Nodal function 6*/
     1334        dbasis[NUMNODESMINI*0+5]=0.;
     1335        dbasis[NUMNODESMINI*1+5]=SQRT3/3.0*(1.0+zeta)/2.0;
     1336        dbasis[NUMNODESMINI*2+5]=0.5*gauss->coord3;
     1337        /*Nodal function 7*/
     1338        dbasis[NUMNODESMINI*0+6]=27.*(1.+zeta)*(1.-zeta)*(-.5*gauss->coord2*gauss->coord3 + .5*gauss->coord1*gauss->coord3);
     1339        dbasis[NUMNODESMINI*1+6]=27.*(1.+zeta)*(1.-zeta)*SQRT3*(-1./6.*gauss->coord2*gauss->coord3 - 1./6.*gauss->coord1*gauss->coord3 +1./3.*gauss->coord1*gauss->coord2);
     1340        dbasis[NUMNODESMINI*2+6]=27*gauss->coord1*gauss->coord2*gauss->coord3*(-2.0*zeta);
    13531341}
    13541342/*}}}*/
    … …  
    13961384/*}}}*/
    13971385/*FUNCTION PentaRef::GetNodalFunctionsP1DerivativesReference {{{*/
    1398 void PentaRef::GetNodalFunctionsP1DerivativesReference(IssmDouble* dl1dl6,GaussPenta* gauss){
     1386void PentaRef::GetNodalFunctionsP1DerivativesReference(IssmDouble* dbasis,GaussPenta* gauss){
    13991387
    14001388        /*This routine returns the values of the nodal functions derivatives  (with respect to the
    14011389         * natural coordinate system) at the gaussian point. Those values vary along xi,eta,z */
    14021390
    1403         IssmDouble A1,A2,A3,z;
    1404 
    1405         A1=gauss->coord1; _assert_(A1>=0 && A1<=1);//first area coordinate value. In term of xi and eta: A1=(1-xi)/2-eta/(2*SQRT3);
    1406         A2=gauss->coord2; _assert_(A2>=0 && A2<=1);//second area coordinate value In term of xi and eta: A2=(1+xi)/2-eta/(2*SQRT3);
    1407         A3=gauss->coord3; _assert_(A3>=0 && A3<=1);//third area coordinate value  In term of xi and eta: A3=y/SQRT3;
    1408         z =gauss->coord4; _assert_(z>=-1 &&  z<=1);//fourth vertical coordinate value. Corresponding nodal function: (1-z)/2 and (1+z)/2
    1409 
    1410         /*First nodal function derivatives. The corresponding nodal function is N=A1*(1-z)/2. Its derivatives follow*/
    1411         *(dl1dl6+NUMNODESP1*0+0)=-0.5*(1.0-z)/2.0;
    1412         *(dl1dl6+NUMNODESP1*1+0)=-0.5/SQRT3*(1.0-z)/2.0;
    1413         *(dl1dl6+NUMNODESP1*2+0)=-0.5*A1;
    1414         //dbasis[NUMNODESP1*0+0]=-0.5*(1.0-zeta)/2.0;
    1415         //dbasis[NUMNODESP1*1+0]=-SQRT3/6.0*(1.0-zeta)/2.0;
    1416         //dbasis[NUMNODESP1*2+0]=-0.5*(-0.5*r-SQRT3/6.0*s+ONETHIRD);
    1417 
    1418         /*Second nodal function: The corresponding nodal function is N=A2*(1-z)/2. Its derivatives follow*/
    1419         *(dl1dl6+NUMNODESP1*0+1)=0.5*(1.0-z)/2.0;
    1420         *(dl1dl6+NUMNODESP1*1+1)=-0.5/SQRT3*(1.0-z)/2.0;
    1421         *(dl1dl6+NUMNODESP1*2+1)=-0.5*A2;
    1422 
    1423         /*Third nodal function: The corresponding nodal function is N=A3*(1-z)/2. Its derivatives follow*/
    1424         *(dl1dl6+NUMNODESP1*0+2)=0.0;
    1425         *(dl1dl6+NUMNODESP1*1+2)=1.0/SQRT3*(1.0-z)/2.0;
    1426         *(dl1dl6+NUMNODESP1*2+2)=-0.5*A3;
    1427 
    1428         /*Fourth nodal function: The corresponding nodal function is N=A1*(1+z)/2. Its derivatives follow*/
    1429         *(dl1dl6+NUMNODESP1*0+3)=-0.5*(1.0+z)/2.0;
    1430         *(dl1dl6+NUMNODESP1*1+3)=-0.5/SQRT3*(1.0+z)/2.0;
    1431         *(dl1dl6+NUMNODESP1*2+3)=0.5*A1;
    1432 
    1433         /*Fifth nodal function: The corresponding nodal function is N=A2*(1+z)/2. Its derivatives follow*/
    1434         *(dl1dl6+NUMNODESP1*0+4)=0.5*(1.0+z)/2.0;
    1435         *(dl1dl6+NUMNODESP1*1+4)=-0.5/SQRT3*(1.0+z)/2.0;
    1436         *(dl1dl6+NUMNODESP1*2+4)=0.5*A2;
    1437 
    1438         /*Sixth nodal function: The corresponding nodal function is N=A3*(1+z)/2. Its derivatives follow*/
    1439         *(dl1dl6+NUMNODESP1*0+5)=0.0;
    1440         *(dl1dl6+NUMNODESP1*1+5)=1.0/SQRT3*(1.0+z)/2.0;
    1441         *(dl1dl6+NUMNODESP1*2+5)=0.5*A3;
     1391        IssmDouble zeta=gauss->coord4;
     1392
     1393        /*Nodal function 1*/
     1394        dbasis[NUMNODESP1*0+0]=-0.5*(1.0-zeta)/2.0;
     1395        dbasis[NUMNODESP1*1+0]=-SQRT3/6.0*(1.0-zeta)/2.0;
     1396        dbasis[NUMNODESP1*2+0]=-0.5*gauss->coord1;
     1397        /*Nodal function 2*/
     1398        dbasis[NUMNODESP1*0+1]=0.5*(1.0-zeta)/2.0;
     1399        dbasis[NUMNODESP1*1+1]=-SQRT3/6.0*(1.0-zeta)/2.0;
     1400        dbasis[NUMNODESP1*2+1]=-0.5*gauss->coord2;
     1401        /*Nodal function 3*/
     1402        dbasis[NUMNODESP1*0+2]=0.;
     1403        dbasis[NUMNODESP1*1+2]=SQRT3/3.0*(1.0-zeta)/2.0;
     1404        dbasis[NUMNODESP1*2+2]=-0.5*gauss->coord3;
     1405        /*Nodal function 4*/
     1406        dbasis[NUMNODESP1*0+3]=-0.5*(1.0+zeta)/2.0;
     1407        dbasis[NUMNODESP1*1+3]=-SQRT3/6.0*(1.0+zeta)/2.0;
     1408        dbasis[NUMNODESP1*2+3]=0.5*gauss->coord1;
     1409        /*Nodal function 5*/
     1410        dbasis[NUMNODESP1*0+4]=0.5*(1.0+zeta)/2.0;
     1411        dbasis[NUMNODESP1*1+4]=-SQRT3/6.0*(1.0+zeta)/2.0;
     1412        dbasis[NUMNODESP1*2+4]=0.5*gauss->coord2;
     1413        /*Nodal function 6*/
     1414        dbasis[NUMNODESP1*0+5]=0.;
     1415        dbasis[NUMNODESP1*1+5]=SQRT3/3.0*(1.0+zeta)/2.0;
     1416        dbasis[NUMNODESP1*2+5]=0.5*gauss->coord3;
    14421417}
    14431418/*}}}*/
Note: See TracChangeset for help on using the changeset viewer.