AcouSTO  version 2.0

◆ build_gmsh()

int build_gmsh ( struct geometry *  geom)

Build nodes from .msh file. Only version >= 2.2 is supported. Only one symmetry plane is supported. Triangular and quadrilateral panels are accepted.

Parameters
[in,out]geomstructure
* Nodes numbering translation scheme.
*
*   v                                        v
*    ^                                      ^
*    |                                      |
*    2                                  1   |
*    |`\                                |`\ |
*    |  `\                              |  `\
*    |    `\                            |   |`\
*    |      `\                          |   +--`\---> u
*    |        `\                        |        `\
*    0----------1--> u                (2=3)--------0
*        GMSH                              ACOUSTO
*
*           v                                v
*          ^                                ^
*          |                                |
*    3-----------2                    2-----------1
*    |     |     |                    |     |     |
*    |     |     |                    |     |     |
*    |     +---- | --> u              |     +-----|---> u
*    |           |                    |           |
*    |           |                    |           |
*    0-----------1                    3-----------0
*        GMSH                            ACOUSTO
*
*
88  {
89  int j, ifreq, inodb,ntmp,netmp,nelmu,nvertu, nome;
90  float ver;
91  long double x,y,z;
92  long double xdotk,dtheta;
93  int i,j1,j2,j3,j4,ftype,datasize,eltype;
94  int ielg,ntag,nphy,ngeo,npart;
95  int icount=0;
96  int ivbase=0;
97  char line[100];
98  FILE* fnodes;
99  struct gmsh_geom *gmshgeom;
100  struct vector xmx0,xsym,proj;
101 
102  char flambda[MAX_PATH];
103  char fgamma[MAX_PATH];
104  char ffunc[MAX_PATH];
105 
106  int indice;
107  struct gmsh_bc* gmshbc;
108 // if (2==modsolinfo->custombc){
109 // file_lambda = fopen(fileinfo->bc_lambda_file,"a");
110 // file_gamma = fopen(fileinfo->bc_gamma_file, "a");
111 // file_f = fopen(fileinfo->bc_func_file, "a");
112 // } // gp end
113 
114  dtheta = 2.0 * acos(-1.0) / runinfo->nsymm;
115  nome = modsolinfo->nome;
116 
117  gmshgeom = (struct gmsh_geom*) geom->ptrgeom;
118  logger(LOG_DEBUG,"loading file %s\n",gmshgeom->filename);
119  /* opening file */
120  fnodes = fopen(gmshgeom->filename,"r");
121  if(NULL == fnodes){
122  logger(LOG_ERROR,MSG_GEOM_FILE_NOTFOUND,gmshgeom->filename);
123  return ACOUSTO_ERROR_GEOMETRY;
124  }
125 
126  /* Reading header */
127  fgets(line,100,fnodes);
128  if(0 != strcmp("$MeshFormat\n",line)){
129  logger(LOG_ERROR,"%s is NOT an ASCII Gmsh file \n",gmshgeom->filename);
130  return ACOUSTO_ERROR_GEOMETRY;
131  }else{
132  logger(LOG_INFO,"%s appears to be an ASCII GMSH file \n",gmshgeom->filename);
133  }
134  fgets(line,100,fnodes);
135  sscanf(line,GMSH_VER_SCANF,&ver,&ftype,&datasize);
136  if(ver < 2.2){
137  logger(LOG_ERROR,"Gmsh file version %f not supported\n",ver);
138  return ACOUSTO_ERROR_GEOMETRY;
139  }
140  if(0 != ftype){
141  logger(LOG_ERROR,"Gmsh file type %d not supported\n",ftype);
142  return ACOUSTO_ERROR_GEOMETRY;
143  }
144 
145 // FILE *file_lambda[nome-1]; // gp start
146 // FILE *file_gamma[nome-1];
147 // FILE *file_f[nome-1];
148  FILE** file_lambda=malloc(sizeof(FILE*)*(nome-1)); // gp start
149  FILE** file_gamma=malloc(sizeof(FILE*)*(nome-1));
150  FILE** file_f=malloc(sizeof(FILE*)*(nome-1));
151 
152  if (2==modsolinfo->custombc) {
153  for (ifreq=0;ifreq<nome;ifreq++) { // B.C. files init
154 
155  sprintf(flambda,fileinfo->bc_lambda_file,ifreq);
156  sprintf(fgamma ,fileinfo->bc_gamma_file ,ifreq);
157  sprintf(ffunc ,fileinfo->bc_func_file ,ifreq);
158 
159  file_lambda[ifreq] = fopen(flambda,"w");
160  file_gamma[ifreq] = fopen(fgamma, "w");
161  file_f[ifreq] = fopen(ffunc, "w");
162  }
163  }
164 
165  icount=0;
166  nelmu = geom->nelmb/runinfo->nsymm;
167  logger(LOG_INFO,"nelmb = %d \n",geom->nelmb);
168  logger(LOG_INFO,"nsymm = %d \n",runinfo->nsymm);
169  logger(LOG_INFO,"nelmu = %d \n",nelmu);
170  nvertu = nelmu*4;
171  while (fgets(line,100,fnodes)!=NULL){
172  if(0 == strcmp("$Nodes\n",line)){
173  fgets(line,100,fnodes);
174  sscanf(line,"%d\n",&ntmp);
175  if(runinfo->nsymm*ntmp != geom->nnodb){
176  logger(LOG_ERROR,"NUMBER OF NODES NOT CONSISTENT %d CHECK MSH FILE\n",ntmp);
177  return ACOUSTO_ERROR_GEOMETRY;
178  }
179  for(inodb=0;inodb<ntmp;inodb++){
180  fgets(line,100,fnodes);
181  sscanf(line,GMSH_NOD_SCANF,&i,&x,&y,&z);
182  geometry->nodes[inodb + _g_inodb].x = x;
183  geometry->nodes[inodb + _g_inodb].y = y;
184  geometry->nodes[inodb + _g_inodb].z = z;
185  }
186  if(2 == runinfo->nsymm){ // Generates symmetrized mesh - 1 symmetry plane
187  logger(LOG_INFO,"GENERATING SYMMETRIZED NODES FOR nsymm = %d \n",runinfo->nsymm);
188  for(inodb=ntmp;inodb<geom->nnodb;inodb++){
189  vec_diff(&xmx0,geometry->nodes[inodb + _g_inodb - ntmp],gmshgeom->refpoint);
190  xdotk=vec_dot(xmx0,gmshgeom->symvec);
191  proj=gmshgeom->symvec;
192  vec_scale(&proj,2.0*xdotk);
193  vec_diff(&xsym,geometry->nodes[inodb + _g_inodb - ntmp],proj);
194  geometry->nodes[inodb + _g_inodb].x = xsym.x;
195  geometry->nodes[inodb + _g_inodb].y = xsym.y;
196  geometry->nodes[inodb + _g_inodb].z = xsym.z;
197  }
198  }
199  if(3 < runinfo->nsymm){ // Generates symmetrized mesh - Axial Symmetry
200  logger(LOG_INFO,"GENERATING SYMMETRIZED NODES FOR nsymm = %d\n",runinfo->nsymm);
201  for(inodb=ntmp;inodb<geom->nnodb;inodb++){
202  vec_copy(&xsym,geometry->nodes[inodb + _g_inodb - ntmp]);
203  vec_copy(&xmx0,gmshgeom->refpoint);
204  vec_copy(&proj,gmshgeom->symvec);
205  vec_rotate2(&xsym, &proj, &xmx0, dtheta);
206  geometry->nodes[inodb + _g_inodb].x = xsym.x;
207  geometry->nodes[inodb + _g_inodb].y = xsym.y;
208  geometry->nodes[inodb + _g_inodb].z = xsym.z;
209  }
210  }
211 
212  }
213  logger(LOG_DEBUG,"LOOP ON NODES FINISHED %d\n",geom->nnodb);
214  if(0 == strcmp("$Elements\n",line)){
215  fgets(line,100,fnodes);
216  sscanf(line,"%d\n",&netmp);
217  for(j=0;j<netmp;j++){
218  fgets(line,100,fnodes);
219  sscanf(line,"%d %d %d\n",&ielg,&eltype,&ntag);
220  if(2==eltype){ // Triangle
221  logger(LOG_DEBUG," Element %d/%d is a TRIANGLE %d\n",j,netmp,eltype);
222  ivbase=icount*4;
223  if(2==ntag){
224  sscanf(line,"%d %d %d %d %d %d %d %d\n",&ielg,&eltype,&ntag,&nphy,&ngeo,&j1,&j2,&j3);
225  geometry->jnodb[ivbase + _g_ielmb] = j2-1;
226  geometry->jnodb[ivbase+1 + _g_ielmb] = j3-1;
227  geometry->jnodb[ivbase+2 + _g_ielmb] = j1-1;
228  geometry->jnodb[ivbase+3 + _g_ielmb] = j1-1;
229  if(2==modsolinfo->custombc){ // gp start
230  logger(LOG_DEBUG," Nphy %d \n",nphy);
231  // trovatag(nphy,&indice);
232  if(-1==trovatag(nphy,&indice)){
233  logger(LOG_ERROR," Physical tags do not match config file. Geometry build fails.\n");
234  return -1;
235  }
236  logger(LOG_DEBUG," Indice %d \n",indice);
237  for (ifreq=0;ifreq<nome;ifreq++) {
238  gmshbc = (struct gmsh_bc*) phystags[ifreq + indice].ptrgmshbc;
239  logger(LOG_DEBUG," Lambda ptr %p gamma ptr %p func ptr %p\n",file_lambda[ifreq],file_gamma[ifreq],file_f[ifreq]);
240  fprintf(file_lambda[ifreq],"%d %f %f\n", ielg, gmshbc->re_lambda, gmshbc->im_lambda);
241  fprintf(file_gamma[ifreq],"%d %f %f\n", ielg, gmshbc->re_gamma, gmshbc->im_gamma);
242  fprintf(file_f[ifreq],"%d %f %f\n", ielg, gmshbc->re_f, gmshbc->im_f);
243  }
244  }
245 // gmshbc = (struct gmsh_bc*) phystags[indice].ptrgmshbc;
246 // fprintf(file_lambda,"%d %f %f\n", ielg, gmshbc->re_lambda, gmshbc->im_lambda);
247 // fprintf(file_gamma,"%d %f %f\n", ielg, gmshbc->re_gamma, gmshbc->im_gamma);
248 // fprintf(file_f,"%d %f %f\n", ielg, gmshbc->re_f, gmshbc->im_f); // gp end
249  }
250  if(3==ntag){
251  logger(LOG_DEBUG," Element %d/%d is a QUADRANGLE %d\n",j,netmp,eltype);
252  sscanf(line,"%d %d %d %d %d %d %d %d %d\n",&ielg,&eltype,&ntag,&nphy,&ngeo,&npart,&j1,&j2,&j3);
253  geometry->jnodb[ivbase + _g_ielmb] = j2-1;
254  geometry->jnodb[ivbase+1 + _g_ielmb] = j3-1;
255  geometry->jnodb[ivbase+2 + _g_ielmb] = j1-1;
256  geometry->jnodb[ivbase+3 + _g_ielmb] = j1-1;
257  if(2==modsolinfo->custombc){ // gp start
258  // trovatag(nphy,&indice);
259  if(-1==trovatag(nphy,&indice)){
260  logger(LOG_ERROR," Physical tags do not match config file. Geometry build fails.\n");
261  return -1;
262  }
263  for (ifreq=0;ifreq<nome;ifreq++) {
264  gmshbc = (struct gmsh_bc*) phystags[ifreq + indice].ptrgmshbc;
265  fprintf(file_lambda[ifreq],"%d %f %f\n", ielg, gmshbc->re_lambda, gmshbc->im_lambda);
266  fprintf(file_gamma[ifreq],"%d %f %f\n", ielg, gmshbc->re_gamma, gmshbc->im_gamma);
267  fprintf(file_f[ifreq],"%d %f %f\n", ielg, gmshbc->re_f, gmshbc->im_f);
268  }
269  }
270 // gmshbc = (struct gmsh_bc*) phystags[indice].ptrgmshbc;
271 // fprintf(file_lambda,"%d %f %f\n", ielg, gmshbc->re_lambda, gmshbc->im_lambda);
272 // fprintf(file_gamma,"%d %f %f\n", ielg, gmshbc->re_gamma, gmshbc->im_gamma);
273 // fprintf(file_f,"%d %f %f\n", ielg, gmshbc->re_f, gmshbc->im_f); // gp end
274  }
275  icount++;
276  }
277  if(3==eltype){ // Quadrangle
278  logger(LOG_DEBUG,"Quadrangles %d\n",eltype);
279  ivbase=icount*4;
280  if(2==ntag){
281  sscanf(line,"%d %d %d %d %d %d %d %d %d\n",&ielg,&eltype,&ntag,&nphy,&ngeo,&j1,&j2,&j3,&j4);
282  geometry->jnodb[ivbase + _g_ielmb] = j2-1;
283  geometry->jnodb[ivbase+1 + _g_ielmb] = j3-1;
284  geometry->jnodb[ivbase+2 + _g_ielmb] = j4-1;
285  geometry->jnodb[ivbase+3 + _g_ielmb] = j1-1;
286  if(2==modsolinfo->custombc){ // gp start
287  if(-1==trovatag(nphy,&indice)){
288  logger(LOG_ERROR," Physical tags do not match config file. Geometry build fails.\n");
289  return -1;
290  }
291  for (ifreq=0;ifreq<nome;ifreq++) {
292  gmshbc = (struct gmsh_bc*) phystags[ifreq + indice].ptrgmshbc;
293  fprintf(file_lambda[ifreq],"%d %f %f\n", ielg, gmshbc->re_lambda, gmshbc->im_lambda);
294  fprintf(file_gamma[ifreq],"%d %f %f\n", ielg, gmshbc->re_gamma, gmshbc->im_gamma);
295  fprintf(file_f[ifreq],"%d %f %f\n", ielg, gmshbc->re_f, gmshbc->im_f);
296  }
297  }
298  }
299  if(3==ntag){
300  sscanf(line,"%d %d %d %d %d %d %d %d %d %d\n",&ielg,&eltype,&ntag,&nphy,&ngeo,&npart,&j1,&j2,&j3,&j4);
301  geometry->jnodb[ivbase + _g_ielmb] = j2-1;
302  geometry->jnodb[ivbase+1 + _g_ielmb] = j3-1;
303  geometry->jnodb[ivbase+2 + _g_ielmb] = j4-1;
304  geometry->jnodb[ivbase+3 + _g_ielmb] = j1-1;
305  if(2==modsolinfo->custombc){ // gp start
306  // trovatag(nphy,&indice);
307  if(-1==trovatag(nphy,&indice)){
308  logger(LOG_ERROR," Physical tags do not match config file. Geometry build fails.\n");
309  return -1;
310  }
311  for (ifreq=0;ifreq<nome;ifreq++) {
312  gmshbc = (struct gmsh_bc*) phystags[ifreq + indice].ptrgmshbc;
313  fprintf(file_lambda[ifreq],"%d %f %f\n", ielg, gmshbc->re_lambda, gmshbc->im_lambda);
314  fprintf(file_gamma[ifreq],"%d %f %f\n", ielg, gmshbc->re_gamma, gmshbc->im_gamma);
315  fprintf(file_f[ifreq],"%d %f %f\n", ielg, gmshbc->re_f, gmshbc->im_f);
316  }
317  }
318  }
319  icount++;
320  }
321  }
322  if(icount != nelmu){
323  logger(LOG_ERROR,"NUMBER OF ELEMENTS NOT CONSISTENT %d CHECK MSH FILE\n",icount);
324  return ACOUSTO_ERROR_GEOMETRY;
325  }
326  }
327  }
328 
329  if (2==modsolinfo->custombc){
330  for (ifreq=0;ifreq<nome;ifreq++) { // closes B.C. files
331  fclose(file_lambda[ifreq]);
332  fclose(file_gamma[ifreq]);
333  fclose(file_f[ifreq]);
334  }
335  }
336 
337 
338  if(2 == runinfo->nsymm){
339  logger(LOG_INFO,"Building symmetrized elements topology - Symmetrized panels circulation inverted. \n",ftype);
340  for(j=nelmu+1;j<=geom->nelmb;j++){
341  ivbase=(j-1)*4;
342  // printf(" %d %d %d %d %d %d\n",geom->nelmb,nelmu,nvertu,j,ivbase,ivbase-nvertu);
343  geometry->jnodb[ivbase + _g_ielmb] = geometry->jnodb[ivbase - nvertu + _g_ielmb] + ntmp;
344  geometry->jnodb[ivbase+1 + _g_ielmb] = geometry->jnodb[ivbase+3 - nvertu + _g_ielmb] + ntmp;
345  geometry->jnodb[ivbase+2 + _g_ielmb] = geometry->jnodb[ivbase+2 - nvertu + _g_ielmb] + ntmp;
346  geometry->jnodb[ivbase+3 + _g_ielmb] = geometry->jnodb[ivbase+1 - nvertu + _g_ielmb] + ntmp;
347  }
348  }
349  if(3 < runinfo->nsymm){
350  logger(LOG_INFO,"Building symmetrized elements topology - Symmetrized panels circulation NOT inverted for AXIALSYMMETRY. \n",ftype);
351  for(j=nelmu+1;j<=geom->nelmb;j++){
352  ivbase=(j-1)*4;
353  // printf(" %d %d %d %d %d %d\n",geom->nelmb,nelmu,nvertu,j,ivbase,ivbase-nvertu);
354  geometry->jnodb[ivbase + _g_ielmb] = geometry->jnodb[ivbase - nvertu + _g_ielmb] + ntmp;
355  geometry->jnodb[ivbase+1 + _g_ielmb] = geometry->jnodb[ivbase+1 - nvertu + _g_ielmb] + ntmp;
356  geometry->jnodb[ivbase+2 + _g_ielmb] = geometry->jnodb[ivbase+2 - nvertu + _g_ielmb] + ntmp;
357  geometry->jnodb[ivbase+3 + _g_ielmb] = geometry->jnodb[ivbase+3 - nvertu + _g_ielmb] + ntmp;
358  }
359  }
360 
361 
362  geom_fill(geom); // filling geometry
363 
364  //adding to global vars
365  _g_icntr += geom->ncntr;
366  _g_inodb += geom->nnodb;
367  _g_ielmb += geom->nelmb;
368 
369  if(gmshgeom->filename) logger(LOG_DEBUG,"file %s read\n",gmshgeom->filename);
370  fclose(fnodes);
371 
372  return 1;
373 
374 }
MSG_GEOM_FILE_NOTFOUND
#define MSG_GEOM_FILE_NOTFOUND
Definition: messages.h:38
gmsh_bc::re_gamma
DOUBLE re_gamma
Definition: structs.h:265
runinfo
struct run_info * runinfo
Definition: globals.h:34
file_info::bc_gamma_file
const char * bc_gamma_file
Definition: structs.h:422
geometry::ncntr
int ncntr
Definition: structs.h:175
gmsh_bc
External geometry boundary conditions assigned with physical tags stored in .msh file using GMSH v2...
Definition: structs.h:261
vec_scale
void vec_scale(struct vector *v1, DOUBLE d)
Definition: math.c:170
GMSH_VER_SCANF
#define GMSH_VER_SCANF
Definition: formats.h:69
geometry::ptrgeom
void * ptrgeom
Definition: structs.h:173
vec_copy
void vec_copy(struct vector *vdest, const struct vector vsrc)
Definition: math.c:158
vector
vector struct to hold triplets.
Definition: structs.h:29
gmsh_bc::im_f
DOUBLE im_f
Definition: structs.h:268
gmsh_geom::refpoint
struct vector refpoint
Definition: structs.h:254
phystags
struct phystag * phystags
Definition: globals.h:58
geom_fill
void geom_fill(struct geometry *geom)
Definition: geom_utils.c:249
LOG_INFO
#define LOG_INFO
Definition: logger.h:26
gmsh_geom::filename
const char * filename
Definition: structs.h:252
fileinfo
struct file_info * fileinfo
Definition: globals.h:32
gmsh_bc::im_lambda
DOUBLE im_lambda
Definition: structs.h:264
vec_rotate2
void vec_rotate2(struct vector *x1, struct vector *nu, struct vector *x0, DOUBLE theta)
Definition: math.c:278
LOG_ERROR
#define LOG_ERROR
Definition: logger.h:24
geometry::nelmb
int nelmb
Definition: structs.h:179
ACOUSTO_ERROR_GEOMETRY
#define ACOUSTO_ERROR_GEOMETRY
Definition: constants.h:25
gmsh_geom::symvec
struct vector symvec
Definition: structs.h:253
modsol_info::custombc
int custombc
Definition: structs.h:523
modsolinfo
struct modsol_info * modsolinfo
Definition: globals.h:44
gmsh_bc::re_f
DOUBLE re_f
Definition: structs.h:267
vec_diff
void vec_diff(struct vector *vdest, const struct vector v1, const struct vector v2)
Definition: math.c:52
vec_dot
DOUBLE vec_dot(struct vector v1, struct vector v2)
Definition: math.c:34
LOG_DEBUG
#define LOG_DEBUG
Definition: logger.h:27
file_info::bc_func_file
const char * bc_func_file
Definition: structs.h:423
malloc
#define malloc(size)
Definition: allocation.h:38
gmsh_geom
External geometry stored in .msh file using GMSH v2.0 format.
Definition: structs.h:250
run_info::nsymm
int nsymm
Definition: structs.h:97
GMSH_NOD_SCANF
#define GMSH_NOD_SCANF
Definition: formats.h:70
gmsh_bc::re_lambda
DOUBLE re_lambda
Definition: structs.h:263
geometry::nnodb
int nnodb
Definition: structs.h:177
gmsh_bc::im_gamma
DOUBLE im_gamma
Definition: structs.h:266
modsol_info::nome
int nome
Definition: structs.h:458
_g_inodb
int _g_inodb
Definition: globals.h:68
file_info::bc_lambda_file
const char * bc_lambda_file
Definition: structs.h:421
_g_ielmb
int _g_ielmb
Definition: globals.h:70
geometry
Geometry info structure.
Definition: structs.h:161
MAX_PATH
#define MAX_PATH
Definition: constants.h:29
trovatag
int trovatag(int phytag, int *indice)
Definition: geom_gmsh.c:32
_g_icntr
int _g_icntr
Definition: globals.h:69
logger
void logger(int level, char *msg,...)
Definition: logger.c:56