Index: ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h =================================================================== --- ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 22074) +++ ../trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 22075) @@ -1066,6 +1066,7 @@ P2xP1Enum, P1xP2Enum, P1xP3Enum, + P1xP4Enum, P2xP4Enum, P1P1Enum, P1P1GLSEnum, Index: ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp =================================================================== --- ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 22074) +++ ../trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 22075) @@ -1021,6 +1021,7 @@ case P2xP1Enum : return "P2xP1"; case P1xP2Enum : return "P1xP2"; case P1xP3Enum : return "P1xP3"; + case P1xP4Enum : return "P1xP4"; case P2xP4Enum : return "P2xP4"; case P1P1Enum : return "P1P1"; case P1P1GLSEnum : return "P1P1GLS"; Index: ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp =================================================================== --- ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 22074) +++ ../trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 22075) @@ -1045,6 +1045,7 @@ else if (strcmp(name,"P2xP1")==0) return P2xP1Enum; else if (strcmp(name,"P1xP2")==0) return P1xP2Enum; else if (strcmp(name,"P1xP3")==0) return P1xP3Enum; + else if (strcmp(name,"P1xP4")==0) return P1xP4Enum; else if (strcmp(name,"P2xP4")==0) return P2xP4Enum; else if (strcmp(name,"P1P1")==0) return P1P1Enum; else if (strcmp(name,"P1P1GLS")==0) return P1P1GLSEnum; Index: ../trunk-jpl/src/c/classes/Elements/PentaRef.cpp =================================================================== --- ../trunk-jpl/src/c/classes/Elements/PentaRef.cpp (revision 22074) +++ ../trunk-jpl/src/c/classes/Elements/PentaRef.cpp (revision 22075) @@ -21,6 +21,7 @@ #define NUMNODESP1b 7 #define NUMNODESP1xP2 9 #define NUMNODESP1xP3 12 +#define NUMNODESP1xP4 15 #define NUMNODESP2xP1 12 #define NUMNODESP2 18 #define NUMNODESP2b 19 @@ -81,6 +82,13 @@ indices[1] = 1; indices[2] = 2; break; + case P1xP4Enum: + numindices = 3; + indices = xNew(numindices); + indices[0] = 0; + indices[1] = 1; + indices[2] = 2; + break; case P2Enum: numindices = 6; indices = xNew(numindices); @@ -422,7 +430,7 @@ basis[ 3]=(9.)/(16.)*gauss->coord1*(zeta-1./3.)*(zeta+1./3.)*(zeta+1.); basis[ 4]=(9.)/(16.)*gauss->coord2*(zeta-1./3.)*(zeta+1./3.)*(zeta+1.); basis[ 5]=(9.)/(16.)*gauss->coord3*(zeta-1./3.)*(zeta+1./3.)*(zeta+1.); - /*quarter-sides of quads*/ + /*third-sides of quads*/ basis[ 6]=(27.)/(16.)*gauss->coord1*(zeta-1)*(zeta-1./3.)*(zeta+1.); basis[ 7]=(27.)/(16.)*gauss->coord2*(zeta-1)*(zeta-1./3.)*(zeta+1.); basis[ 8]=(27.)/(16.)*gauss->coord3*(zeta-1)*(zeta-1./3.)*(zeta+1.); @@ -430,6 +438,26 @@ basis[10]=-(27.)/(16.)*gauss->coord2*(zeta-1)*(zeta+1./3.)*(zeta+1.); basis[11]=-(27.)/(16.)*gauss->coord3*(zeta-1)*(zeta+1./3.)*(zeta+1.); return; + case P1xP4Enum : + /*Corner nodes*/ + basis[ 0]=gauss->coord1*(2./3.)*(zeta-1.)*(zeta-0.5 )*(zeta)*(zeta+0.5); + basis[ 1]=gauss->coord2*(2./3.)*(zeta-1.)*(zeta-0.5 )*(zeta)*(zeta+0.5); + basis[ 2]=gauss->coord3*(2./3.)*(zeta-1.)*(zeta-0.5 )*(zeta)*(zeta+0.5); + basis[ 3]=gauss->coord1*(2./3.)*(zeta-0.5)*(zeta)*(zeta+0.5)*(zeta +1.); + basis[ 4]=gauss->coord2*(2./3.)*(zeta-0.5)*(zeta)*(zeta+0.5)*(zeta +1.); + basis[ 5]=gauss->coord3*(2./3.)*(zeta-0.5)*(zeta)*(zeta+0.5)*(zeta +1.); + /*mid-sides of quads (center of vertical edges)*/ + basis[ 6]=gauss->coord1*4.*(zeta-1.)*(zeta-0.5)*(zeta+0.5)*(zeta+1.); + basis[ 7]=gauss->coord2*4.*(zeta-1.)*(zeta-0.5)*(zeta+0.5)*(zeta+1.); + basis[ 8]=gauss->coord3*4.*(zeta-1.)*(zeta-0.5)*(zeta+0.5)*(zeta+1.); + /*quarter-sides of quads (-0.5 and +0.5 of vertical edges)*/ + basis[ 9]=gauss->coord1*(-8./3.)*(zeta-1.0)*(zeta-0.5)*zeta*(zeta+1.); + basis[10]=gauss->coord2*(-8./3.)*(zeta-1.0)*(zeta-0.5)*zeta*(zeta+1.); + basis[11]=gauss->coord3*(-8./3.)*(zeta-1.0)*(zeta-0.5)*zeta*(zeta+1.); + basis[12]=gauss->coord1*(-8./3.)*(zeta-1.0)*zeta*(zeta+0.5)*(zeta+1.); + basis[13]=gauss->coord2*(-8./3.)*(zeta-1.0)*zeta*(zeta+0.5)*(zeta+1.); + basis[14]=gauss->coord3*(-8./3.)*(zeta-1.0)*zeta*(zeta+0.5)*(zeta+1.); + return; default: _error_("Element type "<coord3 *(27./16.)*( 2.* zeta *( zeta + (1./3.)) + ( zeta - 1. )*( zeta + 1. )); return; + case P1xP4Enum : + /*Nodal function 1*/ + dbasis[NUMNODESP1xP4*0+0 ] = -0.5*(2./3.) *(zeta - 1.)*(zeta - 0.5)*(zeta)*(zeta+0.5); + dbasis[NUMNODESP1xP4*1+0 ] = -SQRT3/6.*(2./3.)*(zeta - 1.)*(zeta - 0.5)*(zeta)*(zeta+0.5); + dbasis[NUMNODESP1xP4*2+0 ] = gauss->coord1 * 2./3.*( (2.*zeta-1)*(zeta -0.5)*(zeta +0.5) + 2.* zeta *zeta *(zeta -1.)); + /*Nodal function 2*/ + dbasis[NUMNODESP1xP4*0+1 ] = +0.5*(2./3.) *(zeta - 1.)*(zeta - 0.5)*(zeta)*(zeta+0.5); + dbasis[NUMNODESP1xP4*1+1 ] = -SQRT3/6.*(2./3.)*(zeta - 1.)*(zeta - 0.5)*(zeta)*(zeta+0.5); + dbasis[NUMNODESP1xP4*2+1 ] = gauss->coord2* 2./3.* ((2.*zeta-1.)*(zeta -0.5)*(zeta +0.5) + 2. * zeta *zeta*(zeta -1.)); + /*Nodal function 3*/ + dbasis[NUMNODESP1xP4*0+2 ] = 0. ; + dbasis[NUMNODESP1xP4*1+2 ] = SQRT3/3.*(2./3.)*(zeta -1.)*(zeta-0.5)*(zeta)*(zeta+0.5); + dbasis[NUMNODESP1xP4*2+2 ] = gauss->coord3* 2./3.*( (2.*zeta-1.)*(zeta -0.5)*(zeta +0.5) + 2.*zeta *zeta*(zeta -1.)); + /*Nodal function 4*/ + dbasis[NUMNODESP1xP4*0+3 ] = -0.5 *(2./3.) *(zeta - 0.5)*(zeta)*(zeta+0.5)*(zeta +1. ); + dbasis[NUMNODESP1xP4*1+3 ] = -SQRT3/6.*(2./3.) *(zeta - 0.5)*(zeta)*(zeta+0.5)*(zeta +1. ); + dbasis[NUMNODESP1xP4*2+3 ] = gauss->coord1* 2./3.*( (2.*zeta+1.)*(zeta -0.5)*(zeta +0.5) + 2.*zeta *zeta*(zeta +1.)); + /*Nodal function 5*/ + dbasis[NUMNODESP1xP4*0+4 ] = +0.5 * (2./3.) *(zeta - 0.5)*(zeta)*(zeta+0.5)*(zeta +1. ); + dbasis[NUMNODESP1xP4*1+4 ] = -SQRT3/6.*(2./3.)*(zeta - 0.5)*(zeta)*(zeta+0.5)*(zeta +1. ); + dbasis[NUMNODESP1xP4*2+4 ] = gauss->coord2 * 2./3.*( (2.*zeta+1.)*(zeta -0.5)*(zeta +0.5) + 2.*zeta *zeta*(zeta +1.)); + /*Nodal function 6*/ + dbasis[NUMNODESP1xP4*0+5 ] = 0. ; + dbasis[NUMNODESP1xP4*1+5 ] = SQRT3/3.*(2./3.) *(zeta - 0.5)*(zeta)*(zeta+0.5)*(zeta +1. ); + dbasis[NUMNODESP1xP4*2+5 ] = gauss->coord3 * 2./3.*( (2.*zeta+1.)*(zeta -0.5)*(zeta +0.5) + 2.*zeta *zeta*(zeta +1)); + /*Nodal function 7*/ + dbasis[NUMNODESP1xP4*0+6 ] = -0.5 * 4.*(zeta - 1.)*(zeta - 0.5)*(zeta+0.5)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*1+6 ] = -SQRT3/6.* 4.*(zeta - 1.)*(zeta - 0.5)*(zeta+0.5)*(zeta +1.) ; + dbasis[NUMNODESP1xP4*2+6 ] = gauss->coord1* 4.*( 4.*zeta *zeta*zeta - (5./2.)*zeta ); + /*Nodal function 8*/ + dbasis[NUMNODESP1xP4*0+7 ] = +0.5* 4.*(zeta - 1.)*(zeta - 0.5)*(zeta+0.5)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*1+7 ] = -SQRT3/6.* 4.*(zeta - 1.)*(zeta - 0.5)*(zeta+0.5)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*2+7 ] = gauss->coord2* 4.*( 4.*zeta *zeta*zeta - (5./2.)*zeta ); + + /*Nodal function 9*/ + dbasis[NUMNODESP1xP4*0+8 ] = 0. ; + dbasis[NUMNODESP1xP4*1+8 ] = SQRT3/3. * 4.*(zeta - 1.)*(zeta - 0.5)*(zeta+0.5)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*2+8 ] = gauss->coord3* 4.*( 4.*zeta *zeta*zeta - (5./2.)*zeta ); + + /*Nodal function 10*/ + dbasis[NUMNODESP1xP4*0+9 ] = -0.5* (-8./3.)*(zeta - 1.)*(zeta - 0.5)*(zeta)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*1+9 ] = -SQRT3/3.* (-8./3.)*(zeta - 1.)*(zeta - 0.5)*(zeta)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*2+9 ] = gauss->coord1* (-8./3.)*((2.*zeta -1.)*(zeta-0.5)*(zeta +1.) +zeta*(zeta -1.)*( 2.*zeta + 0.5)); + /*Nodal function 11*/ + dbasis[NUMNODESP1xP4*0+10] = +0.5* (-8./3.)*(zeta - 1.)*(zeta - 0.5)*(zeta)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*1+10] = -SQRT3/6.* (-8./3.)*(zeta - 1.)*(zeta - 0.5)*(zeta)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*2+10] = gauss->coord2* (-8./3.)*((2.*zeta -1.)*(zeta-0.5)*(zeta +1.) +zeta *(zeta -1.)*( 2.*zeta + 0.5)); + /*Nodal function 12*/ + dbasis[NUMNODESP1xP4*0+11] = 0. ; + dbasis[NUMNODESP1xP4*1+11] = SQRT3/6.* (-8./3.)*(zeta - 1.)*(zeta - 0.5)*(zeta)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*2+11] = gauss->coord3* (-8./3.)*((2.*zeta-1.)*(zeta-0.5)*(zeta +1.) +zeta *(zeta -1.)*( 2.*zeta + 0.5)); + /*Nodal function 13*/ + dbasis[NUMNODESP1xP4*0+12] = -0.5* (-8./3.)*(zeta - 1.)*(zeta)*(zeta+0.5)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*1+12] = -SQRT3/6.* (-8./3.)*(zeta - 1.)*(zeta)*(zeta+0.5)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*2+12] = gauss->coord1* (-8./3.)*((2.*zeta -1. ) *(zeta+0.5)* (zeta +1.) + zeta* (zeta -1.)*( 2.*zeta + 3./2.)); + /*Nodal function 14*/ + dbasis[NUMNODESP1xP4*0+13] = +0.5* (-8./3.)*(zeta - 1.)*(zeta)*(zeta+0.5)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*1+13] = -SQRT3/6.* (-8./3.)*(zeta - 1.)*(zeta)*(zeta+0.5)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*2+13] = gauss->coord2* (-8./3.)*((2.*zeta -1. )*(zeta+0.5)*(zeta +1.) + zeta*(zeta -1.)*( 2.*zeta + 3./2.)); + /*Nodal function 15*/ + dbasis[NUMNODESP1xP4*0+14] = 0 ; + dbasis[NUMNODESP1xP4*1+14] = SQRT3/6.* (-8./3.)*(zeta - 1.)*(zeta + 0.5)*(zeta)*(zeta +1. ) ; + dbasis[NUMNODESP1xP4*2+14] = gauss->coord3* (-8./3.)*((2. *zeta -1. )*(zeta+0.5)*(zeta +1.) + zeta*(zeta -1.)*( 2.*zeta + 3./2.)); + + return; default: _error_("Element type "<(numseg*2); + indices[0*2 + 0] = 0; indices[0*2 + 1] = 6; + indices[1*2 + 0] = 1; indices[1*2 + 1] = 7; + indices[2*2 + 0] = 2; indices[2*2 + 1] = 8; + break; case P2Enum: numseg = 6; indices = xNew(numseg*2); @@ -1129,6 +1229,7 @@ case P1xP2Enum: return NUMNODESP1xP2; case P2xP4Enum: return NUMNODESP2xP4; case P1xP3Enum: return NUMNODESP1xP3; + case P1xP4Enum: return NUMNODESP1xP4; case P1P1Enum: return NUMNODESP1*2; case P1P1GLSEnum: return NUMNODESP1*2; case MINIcondensedEnum: return NUMNODESP1b+NUMNODESP1; @@ -1195,6 +1296,7 @@ break; case P1xP2Enum: case P1xP3Enum: + case P1xP4Enum: numindices = 3; indices = xNew(numindices); indices[0] = 3; Index: ../trunk-jpl/src/c/classes/gauss/GaussPenta.cpp =================================================================== --- ../trunk-jpl/src/c/classes/gauss/GaussPenta.cpp (revision 22074) +++ ../trunk-jpl/src/c/classes/gauss/GaussPenta.cpp (revision 22075) @@ -600,6 +600,29 @@ default: _error_("node index should be in [0 11]"); } break; + case P1xP4Enum: + switch(iv){ + case 0 : coord1=1.; coord2=0.; coord3=0.; coord4=-1.; break; + case 1 : coord1=0.; coord2=1.; coord3=0.; coord4=-1.; break; + case 2 : coord1=0.; coord2=0.; coord3=1.; coord4=-1.; break; + case 3 : coord1=1.; coord2=0.; coord3=0.; coord4=+1.; break; + case 4 : coord1=0.; coord2=1.; coord3=0.; coord4=+1.; break; + case 5 : coord1=0.; coord2=0.; coord3=1.; coord4=+1.; break; + + case 6 : coord1=1.; coord2=0.; coord3=0.; coord4=0.; break; + case 7 : coord1=0.; coord2=1.; coord3=0.; coord4=0.; break; + case 8 : coord1=0.; coord2=0.; coord3=1.; coord4=0.; break; + + case 9 : coord1=1.; coord2=0.; coord3=0.; coord4=-0.5; break; + case 10: coord1=0.; coord2=1.; coord3=0.; coord4=-0.5; break; + case 11: coord1=0.; coord2=0.; coord3=1.; coord4=-0.5; break; + + case 12: coord1=1.; coord2=0.; coord3=0.; coord4=+0.5; break; + case 13: coord1=0.; coord2=1.; coord3=0.; coord4=+0.5; break; + case 14: coord1=0.; coord2=0.; coord3=1.; coord4=+0.5; break; + default: _error_("node index should be in [0 14]"); + } + break; case P2xP1Enum: switch(iv){ case 0: coord1=1.; coord2=0.; coord3=0.; coord4=-1.; break; Index: ../trunk-jpl/src/c/classes/Inputs/PentaInput.cpp =================================================================== --- ../trunk-jpl/src/c/classes/Inputs/PentaInput.cpp (revision 22074) +++ ../trunk-jpl/src/c/classes/Inputs/PentaInput.cpp (revision 22075) @@ -253,6 +253,18 @@ for(int i=2;i<9;i+=3) this->values[i]=this->values[11]; } break; + case P1xP4Enum: + if(start==-1){ + for(int i=0+3;i<15;i+=3) this->values[i]=this->values[0]; + for(int i=1+3;i<15;i+=3) this->values[i]=this->values[1]; + for(int i=2+3;i<15;i+=3) this->values[i]=this->values[2]; + } + else{ + for(int i=0;i<12;i+=3) this->values[i]=this->values[12]; + for(int i=1;i<12;i+=3) this->values[i]=this->values[13]; + for(int i=2;i<12;i+=3) this->values[i]=this->values[14]; + } + break; default: _error_("not supported yet for type "<interpolation_type)); }