source: issm/trunk-jpl/src/c/shared/Elements/CoordinateSystemTransform.cpp@ 12493

Last change on this file since 12493 was 12493, checked in by Mathieu Morlighem, 13 years ago

replaced all _error_ to _error2_
File size: 2.8 KB
Line 
1/*!\file: GetGlobalDofList.cpp
2 * \brief create transform matrix for different coordinate systems
3 */
4#include "./elements.h"
5#include <math.h>
6
7void CoordinateSystemTransform(IssmDouble** ptransform,Node** nodes,int numnodes,int* cs_array){
8
9 int i,counter;
10 int numdofs = 0;
11 IssmDouble norm;
12 IssmDouble *transform = NULL;
13 IssmDouble *values = NULL;
14 IssmDouble coord_system[3][3];
15
16 /*Some checks in debugging mode*/
17 _assert_(numnodes && nodes);
18
19 /*Get total number of dofs*/
20 for(i=0;i<numnodes;i++){
21 switch(cs_array[i]){
22 case XYEnum: numdofs+=2; break;
23 case XYZPEnum: numdofs+=4; break;
24 default: _error2_("Coordinate system " << EnumToStringx(cs_array[i]) << " not supported yet");
25 }
26 }
27
28 /*Allocate and initialize transform matrix*/
29 transform=xNew<IssmDouble>(numdofs*numdofs);
30 for(i=0;i<numdofs*numdofs;i++) transform[i]=0.0;
31
32 /*Create transform matrix for all nodes (x,y for 2d and x,y,z for 3d). It is a block matrix
33 *for 3 nodes:
34
35 * | T1 0 0 |
36 * Q = | 0 T2 0 |
37 * | 0 0 T3|
38 *
39 * Where T1 is the transform matrix for node 1. It is a simple copy of the coordinate system
40 * associated to this node*/
41 counter=0;
42 for(i=0;i<numnodes;i++){
43 nodes[i]->GetCoordinateSystem(&coord_system[0][0]);
44 switch(cs_array[i]){
45 case XYEnum:
46 /*We remove the z component, we need to renormalize x and y: x=[x1 x2 0] y=[-x2 x1 0]*/
47 norm = sqrt( coord_system[0][0]*coord_system[0][0] + coord_system[1][0]*coord_system[1][0]); _assert_(norm>1.e-4);
48 transform[(numdofs)*(counter+0) + counter+0] = coord_system[0][0]/norm;
49 transform[(numdofs)*(counter+0) + counter+1] = - coord_system[1][0]/norm;
50 transform[(numdofs)*(counter+1) + counter+0] = coord_system[1][0]/norm;
51 transform[(numdofs)*(counter+1) + counter+1] = coord_system[0][0]/norm;
52 counter+=2;
53 break;
54 case XYZPEnum:
55 /*Only the first 3 coordinates are changed (x,y,z), leave the others (P) unchanged*/
56 transform[(numdofs)*(counter+0) + counter+0] = coord_system[0][0];
57 transform[(numdofs)*(counter+0) + counter+1] = coord_system[0][1];
58 transform[(numdofs)*(counter+0) + counter+2] = coord_system[0][2];
59 transform[(numdofs)*(counter+1) + counter+0] = coord_system[1][0];
60 transform[(numdofs)*(counter+1) + counter+1] = coord_system[1][1];
61 transform[(numdofs)*(counter+1) + counter+2] = coord_system[1][2];
62 transform[(numdofs)*(counter+2) + counter+0] = coord_system[2][0];
63 transform[(numdofs)*(counter+2) + counter+1] = coord_system[2][1];
64 transform[(numdofs)*(counter+2) + counter+2] = coord_system[2][2];
65 transform[(numdofs)*(counter+3) + counter+3] = 1.0;
66 counter+=4;
67 break;
68 default:
69 _error2_("Coordinate system " << EnumToStringx(cs_array[i]) << " not supported yet");
70 }
71 }
72
73 /*Assign output pointer*/
74 *ptransform=transform;
75}
Note: See TracBrowser for help on using the repository browser.