Changeset 4769 for issm/trunk/src/c/objects/Inputs/TriaVertexInput.cpp
- Timestamp:
- 07/22/10 15:46:39 (15 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
issm/trunk/src/c/objects/Inputs/TriaVertexInput.cpp
r4546 r4769 24 24 /*}}}*/ 25 25 /*FUNCTION TriaVertexInput::TriaVertexInput(int in_enum_type,double* values){{{1*/ 26 TriaVertexInput::TriaVertexInput(int in_enum_type,double* in_values){ 27 26 TriaVertexInput::TriaVertexInput(int in_enum_type,double* in_values) 27 :TriaRef(1) 28 { 29 30 /*Set TriaRef*/ 31 this->SetElementType(P1Enum,0); 32 this->element_type=P1Enum; 33 34 /*Set Enum*/ 28 35 enum_type=in_enum_type; 36 37 /*Set values*/ 29 38 values[0]=in_values[0]; 30 39 values[1]=in_values[1]; … … 311 320 312 321 /*Intermediary*/ 313 /*FUNCTION TriaVertexInput::GetNodalFunctions {{{1*/314 void TriaVertexInput::GetNodalFunctions(double* l1l2l3, double* gauss){315 316 /*This routine returns the values of the nodal functions at the gaussian point.*/317 318 /*First nodal function: */319 l1l2l3[0]=gauss[0];320 321 /*Second nodal function: */322 l1l2l3[1]=gauss[1];323 324 /*Third nodal function: */325 l1l2l3[2]=gauss[2];326 327 }328 /*}}}*/329 /*FUNCTION TriaVertexInput::GetB {{{1*/330 void TriaVertexInput::GetB(double* B, double* xyz_list, double* gauss){331 /*Compute B matrix. B=[B1 B2 B3] where Bi is of size 3*NDOF2.332 * For grid i, Bi can be expressed in the actual coordinate system333 * by:334 * Bi=[ dh/dx 0 ]335 * [ 0 dh/dy ]336 * [ 1/2*dh/dy 1/2*dh/dx ]337 * where h is the interpolation function for grid i.338 *339 * We assume B has been allocated already, of size: 3x(NDOF2*numgrids)340 */341 342 int i;343 const int NDOF2=2;344 const int numgrids=3;345 346 double dh1dh3[NDOF2][numgrids];347 348 349 /*Get dh1dh2dh3 in actual coordinate system: */350 GetNodalFunctionsDerivatives(&dh1dh3[0][0],xyz_list,gauss);351 352 /*Build B: */353 for (i=0;i<numgrids;i++){354 *(B+NDOF2*numgrids*0+NDOF2*i)=dh1dh3[0][i]; //B[0][NDOF2*i]=dh1dh3[0][i];355 *(B+NDOF2*numgrids*0+NDOF2*i+1)=0;356 *(B+NDOF2*numgrids*1+NDOF2*i)=0;357 *(B+NDOF2*numgrids*1+NDOF2*i+1)=dh1dh3[1][i];358 *(B+NDOF2*numgrids*2+NDOF2*i)=(float).5*dh1dh3[1][i];359 *(B+NDOF2*numgrids*2+NDOF2*i+1)=(float).5*dh1dh3[0][i];360 }361 }362 /*}}}*/363 /*FUNCTION TriaVertexInput::GetNodalFunctionsDerivatives {{{1*/364 void TriaVertexInput::GetNodalFunctionsDerivatives(double* dh1dh3,double* xyz_list, double* gauss){365 366 /*This routine returns the values of the nodal functions derivatives (with respect to the367 * actual coordinate system: */368 int i;369 const int NDOF2=2;370 const int numgrids=3;371 double dh1dh3_ref[NDOF2][numgrids];372 double Jinv[NDOF2][NDOF2];373 374 /*Get derivative values with respect to parametric coordinate system: */375 GetNodalFunctionsDerivativesReference(&dh1dh3_ref[0][0], gauss);376 377 /*Get Jacobian invert: */378 GetJacobianInvert(&Jinv[0][0], xyz_list, gauss);379 380 /*Build dh1dh3:381 *382 * [dhi/dx]= Jinv*[dhi/dr]383 * [dhi/dy] [dhi/ds]384 */385 for (i=0;i<numgrids;i++){386 dh1dh3[numgrids*0+i]=Jinv[0][0]*dh1dh3_ref[0][i]+Jinv[0][1]*dh1dh3_ref[1][i];387 dh1dh3[numgrids*1+i]=Jinv[1][0]*dh1dh3_ref[0][i]+Jinv[1][1]*dh1dh3_ref[1][i];388 }389 390 }391 /*}}}*/392 /*FUNCTION TriaVertexInput::GetNodalFunctionsDerivativesReference {{{1*/393 void TriaVertexInput::GetNodalFunctionsDerivativesReference(double* dl1dl3,double* gauss_l1l2l3){394 395 /*This routine returns the values of the nodal functions derivatives (with respect to the396 * natural coordinate system) at the gaussian point. */397 398 const int NDOF2=2;399 const int numgrids=3;400 401 /*First nodal function: */402 *(dl1dl3+numgrids*0+0)=-0.5;403 *(dl1dl3+numgrids*1+0)=-1.0/(2.0*SQRT3);404 405 /*Second nodal function: */406 *(dl1dl3+numgrids*0+1)=0.5;407 *(dl1dl3+numgrids*1+1)=-1.0/(2.0*SQRT3);408 409 /*Third nodal function: */410 *(dl1dl3+numgrids*0+2)=0;411 *(dl1dl3+numgrids*1+2)=1.0/SQRT3;412 413 }414 /*}}}*/415 /*FUNCTION TriaVertexInput::GetJacobian {{{1*/416 void TriaVertexInput::GetJacobian(double* J, double* xyz_list,double* gauss){417 418 /*The Jacobian is constant over the element, discard the gaussian points.419 * J is assumed to have been allocated of size NDOF2xNDOF2.*/420 421 const int NDOF2=2;422 const int numgrids=3;423 double x1,y1,x2,y2,x3,y3;424 425 x1=*(xyz_list+numgrids*0+0);426 y1=*(xyz_list+numgrids*0+1);427 x2=*(xyz_list+numgrids*1+0);428 y2=*(xyz_list+numgrids*1+1);429 x3=*(xyz_list+numgrids*2+0);430 y3=*(xyz_list+numgrids*2+1);431 432 433 *(J+NDOF2*0+0)=0.5*(x2-x1);434 *(J+NDOF2*1+0)=SQRT3/6.0*(2*x3-x1-x2);435 *(J+NDOF2*0+1)=0.5*(y2-y1);436 *(J+NDOF2*1+1)=SQRT3/6.0*(2*y3-y1-y2);437 }438 /*}}}*/439 /*FUNCTION TriaVertexInput::GetJacobianInvert {{{1*/440 void TriaVertexInput::GetJacobianInvert(double* Jinv, double* xyz_list,double* gauss){441 442 double Jdet;443 const int NDOF2=2;444 const int numgrids=3;445 446 /*Call Jacobian routine to get the jacobian:*/447 GetJacobian(Jinv, xyz_list, gauss);448 449 /*Invert Jacobian matrix: */450 MatrixInverse(Jinv,NDOF2,NDOF2,NULL,0,&Jdet);451 452 }453 /*}}}*/454 322 /*FUNCTION TriaVertexInput::SquareMin{{{1*/ 455 323 void TriaVertexInput::SquareMin(double* psquaremin, bool process_units,Parameters* parameters){
Note:
See TracChangeset
for help on using the changeset viewer.