CalculiX  2.8
A Free Software Three-Dimensional Structural Finite Element Program
 All Classes Files Functions Variables Macros
mastructcs.c File Reference
#include <stdlib.h>
#include <math.h>
#include <stdio.h>
#include <string.h>
#include "CalculiX.h"
Include dependency graph for mastructcs.c:

Go to the source code of this file.

Macros

#define min(a, b)   ((a) <= (b) ? (a) : (b))
 
#define max(a, b)   ((a) >= (b) ? (a) : (b))
 

Functions

void mastructcs (ITG *nk, ITG *kon, ITG *ipkon, char *lakon, ITG *ne, ITG *nodeboun, ITG *ndirboun, ITG *nboun, ITG *ipompc, ITG *nodempc, ITG *nmpc, ITG *nactdof, ITG *icol, ITG *jq, ITG **mast1p, ITG **irowp, ITG *isolver, ITG *neq, ITG *ikmpc, ITG *ilmpc, ITG *ipointer, ITG *nzs, ITG *nmethod, ITG *ics, double *cs, char *labmpc, ITG *mcs, ITG *mi, ITG *mortar)
 

Macro Definition Documentation

#define max (   a,
 
)    ((a) >= (b) ? (a) : (b))

Definition at line 25 of file mastructcs.c.

#define min (   a,
 
)    ((a) <= (b) ? (a) : (b))

Definition at line 24 of file mastructcs.c.

Function Documentation

void mastructcs ( ITG nk,
ITG kon,
ITG ipkon,
char *  lakon,
ITG ne,
ITG nodeboun,
ITG ndirboun,
ITG nboun,
ITG ipompc,
ITG nodempc,
ITG nmpc,
ITG nactdof,
ITG icol,
ITG jq,
ITG **  mast1p,
ITG **  irowp,
ITG isolver,
ITG neq,
ITG ikmpc,
ITG ilmpc,
ITG ipointer,
ITG nzs,
ITG nmethod,
ITG ics,
double *  cs,
char *  labmpc,
ITG mcs,
ITG mi,
ITG mortar 
)

Definition at line 27 of file mastructcs.c.

33  {
34 
35  /* determines the structure of the thermo-mechanical matrices with
36  cyclic symmetry;
37  (i.e. the location of the nonzeros */
38 
39  char lakonl[2]=" \0";
40 
41  ITG i,j,k,l,jj,ll,id,index,jdof1,jdof2,idof1,idof2,mpc1,mpc2,id1,id2,
42  ist1,ist2,node1,node2,isubtract,nmast,ifree,istart,istartold,
43  index1,index2,m,node,nzs_,ist,kflag,indexe,nope,isize,*mast1=NULL,
44  *irow=NULL,inode,icomplex,inode1,icomplex1,inode2,
45  icomplex2,kdof1,kdof2,ilength,lprev,ij,mt=mi[1]+1;
46 
47  /* the indices in the comments follow FORTRAN convention, i.e. the
48  fields start with 1 */
49 
50  mast1=*mast1p;
51  irow=*irowp;
52 
53  kflag=2;
54  nzs_=nzs[1];
55 
56  /* initialisation of nactmpc */
57 
58  for(i=0;i<mt**nk;++i){nactdof[i]=0;}
59 
60  /* determining the active degrees of freedom due to elements */
61 
62  for(i=0;i<*ne;++i){
63 
64  if(ipkon[i]<0) continue;
65  indexe=ipkon[i];
66 /* Bernhardi start */
67  if (strcmp1(&lakon[8*i+3],"8I")==0)nope=11;
68  else if(strcmp1(&lakon[8*i+3],"20")==0)nope=20;
69 /* Bernhardi end */
70  else if(strcmp1(&lakon[8*i+3],"2")==0)nope=26;
71  else if (strcmp1(&lakon[8*i+3],"8")==0)nope=8;
72  else if (strcmp1(&lakon[8*i+3],"10")==0)nope=10;
73  else if ((strcmp1(&lakon[8*i+3],"4")==0)||
74  (strcmp1(&lakon[8*i+2],"4")==0)) nope=4;
75  else if (strcmp1(&lakon[8*i+3],"15")==0)nope=15;
76  else if (strcmp1(&lakon[8*i+3],"6")==0)nope=6;
77  else if (strcmp1(&lakon[8*i],"E")==0){
78  if((strcmp1(&lakon[8*i+6],"C")==0)&&(*mortar==1)){
79  nope=kon[ipkon[i]-1];
80  }else{
81  lakonl[0]=lakon[8*i+7];
82  nope=atoi(lakonl)+1;
83  }
84  }else continue;
85 
86 /* else if (strcmp1(&lakon[8*i],"E")==0){
87  lakonl[0]=lakon[8*i+7];
88  nope=atoi(lakonl)+1;}
89  else continue;*/
90 
91  for(j=0;j<nope;++j){
92  node=kon[indexe+j]-1;
93  for(k=1;k<4;++k){
94  nactdof[mt*node+k]=1;
95  }
96  }
97  }
98 
99  /* determining the active degrees of freedom due to mpc's */
100 
101  for(i=0;i<*nmpc;++i){
102  index=ipompc[i]-1;
103  do{
104  if((nodempc[3*index+1]!=0)&&(nodempc[3*index+1]<4)){
105  nactdof[mt*(nodempc[3*index]-1)+nodempc[3*index+1]]=1;}
106  index=nodempc[3*index+2];
107  if(index==0) break;
108  index--;
109  }while(1);
110  }
111 
112  /* subtracting the SPC and MPC nodes */
113 
114  for(i=0;i<*nboun;++i){
115  if(ndirboun[i]>mi[1]) continue;
116  nactdof[mt*(nodeboun[i]-1)+ndirboun[i]]=0;
117  }
118 
119  for(i=0;i<*nmpc;++i){
120  index=ipompc[i]-1;
121  if(nodempc[3*index+1]>mi[1]) continue;
122  nactdof[mt*(nodempc[3*index]-1)+nodempc[3*index+1]]=0;
123  }
124 
125  /* numbering the active degrees of freedom */
126 
127  neq[0]=0;
128  for(i=0;i<*nk;++i){
129  for(j=1;j<4;++j){
130  if(nactdof[mt*i+j]!=0){
131  ++neq[0];
132  nactdof[mt*i+j]=neq[0];
133  }
134  }
135  }
136 
137  ifree=0;
138 
139  /* determining the position of each nonzero matrix element
140 
141  mast1(ipointer(i)) = first nonzero row in column i
142  irow(ipointer(i)) points to further nonzero elements in
143  column i */
144 
145  for(i=0;i<6**nk;++i){ipointer[i]=0;}
146 
147  for(i=0;i<*ne;++i){
148 
149  if(ipkon[i]<0) continue;
150  indexe=ipkon[i];
151 /* Bernhardi start */
152  if(strcmp1(&lakon[8*i],"C3D8I")==0){nope=11;}
153  else if(strcmp1(&lakon[8*i+3],"20")==0)nope=20;
154 /* Bernhardi end */
155  else if(strcmp1(&lakon[8*i+3],"2")==0)nope=26;
156  else if (strcmp1(&lakon[8*i+3],"8")==0)nope=8;
157  else if (strcmp1(&lakon[8*i+3],"10")==0)nope=10;
158  else if (strcmp1(&lakon[8*i+3],"4")==0)nope=4;
159  else if (strcmp1(&lakon[8*i+3],"15")==0)nope=15;
160  else if (strcmp1(&lakon[8*i+3],"6")==0)nope=6;
161  else if (strcmp1(&lakon[8*i],"E")==0){
162  if((strcmp1(&lakon[8*i+6],"C")==0)&&(*mortar==1)){
163  nope=kon[ipkon[i]-1];
164  }else{
165  lakonl[0]=lakon[8*i+7];
166  nope=atoi(lakonl)+1;
167  }
168  }else continue;
169 
170 /* else if (strcmp1(&lakon[8*i],"E")==0){
171  lakonl[0]=lakon[8*i+7];
172  nope=atoi(lakonl)+1;}
173  else continue;*/
174 
175  for(jj=0;jj<3*nope;++jj){
176 
177  j=jj/3;
178  k=jj-3*j;
179 
180  node1=kon[indexe+j];
181  jdof1=nactdof[mt*(node1-1)+k+1];
182 
183  for(ll=jj;ll<3*nope;++ll){
184 
185  l=ll/3;
186  m=ll-3*l;
187 
188  node2=kon[indexe+l];
189  jdof2=nactdof[mt*(node2-1)+m+1];
190 
191  /* check whether one of the DOF belongs to a SPC or MPC */
192 
193  if((jdof1!=0)&&(jdof2!=0)){
194  insert(ipointer,&mast1,&irow,&jdof1,&jdof2,&ifree,&nzs_);
195  kdof1=jdof1+neq[0];kdof2=jdof2+neq[0];
196  insert(ipointer,&mast1,&irow,&kdof1,&kdof2,&ifree,&nzs_);
197  }
198  else if((jdof1!=0)||(jdof2!=0)){
199 
200  /* idof1: genuine DOF
201  idof2: nominal DOF of the SPC/MPC */
202 
203  if(jdof1==0){
204  idof1=jdof2;
205  idof2=8*node1+k-7;}
206  else{
207  idof1=jdof1;
208  idof2=8*node2+m-7;}
209 
210  if(*nmpc>0){
211 
212  FORTRAN(nident,(ikmpc,&idof2,nmpc,&id));
213  if((id>0)&&(ikmpc[id-1]==idof2)){
214 
215  /* regular DOF / MPC */
216 
217  id1=ilmpc[id-1];
218  ist=ipompc[id1-1];
219  index=nodempc[3*ist-1];
220  if(index==0) continue;
221  while(1){
222  inode=nodempc[3*index-3];
223  icomplex=0;
224  if(strcmp1(&labmpc[(id1-1)*20],"CYCLIC")==0){
225  icomplex=atoi(&labmpc[20*(id1-1)+6]);
226  }
227  else if(strcmp1(&labmpc[(id1-1)*20],"SUBCYCLIC")==0){
228  for(ij=0;ij<*mcs;ij++){
229  ilength=cs[17*ij+3];
230  lprev=cs[17*ij+13];
231  FORTRAN(nident,(&ics[lprev],&inode,&ilength,&id));
232  if(id>0){
233  if(ics[lprev+id-1]==inode){
234  icomplex=ij+1;
235  break;
236  }
237  }
238  }
239  }
240 // idof2=nactdof[mt*inode+nodempc[3*index-2]-4];
241  idof2=nactdof[mt*(inode-1)+nodempc[3*index-2]];
242  if(idof2!=0){
243  insert(ipointer,&mast1,&irow,&idof1,&idof2,&ifree,&nzs_);
244  kdof1=idof1+neq[0];kdof2=idof2+neq[0];
245  insert(ipointer,&mast1,&irow,&kdof1,&kdof2,&ifree,&nzs_);
246  if((icomplex!=0)&&(idof1!=idof2)){
247  insert(ipointer,&mast1,&irow,&kdof1,&idof2,&ifree,&nzs_);
248  insert(ipointer,&mast1,&irow,&idof1,&kdof2,&ifree,&nzs_);
249  }
250  }
251  index=nodempc[3*index-1];
252  if(index==0) break;
253  }
254  continue;
255  }
256  }
257  }
258 
259  else{
260  idof1=8*node1+k-7;
261  idof2=8*node2+m-7;
262  mpc1=0;
263  mpc2=0;
264  if(*nmpc>0){
265  FORTRAN(nident,(ikmpc,&idof1,nmpc,&id1));
266  if((id1>0)&&(ikmpc[id1-1]==idof1)) mpc1=1;
267  FORTRAN(nident,(ikmpc,&idof2,nmpc,&id2));
268  if((id2>0)&&(ikmpc[id2-1]==idof2)) mpc2=1;
269  }
270  if((mpc1==1)&&(mpc2==1)){
271  id1=ilmpc[id1-1];
272  id2=ilmpc[id2-1];
273  if(id1==id2){
274 
275  /* MPC id1 / MPC id1 */
276 
277  ist=ipompc[id1-1];
278  index1=nodempc[3*ist-1];
279  if(index1==0) continue;
280  while(1){
281  inode1=nodempc[3*index1-3];
282  icomplex1=0;
283  if(strcmp1(&labmpc[(id1-1)*20],"CYCLIC")==0){
284  icomplex1=atoi(&labmpc[20*(id1-1)+6]);
285  }
286  else if(strcmp1(&labmpc[(id1-1)*20],"SUBCYCLIC")==0){
287  for(ij=0;ij<*mcs;ij++){
288  ilength=cs[17*ij+3];
289  lprev=cs[17*ij+13];
290  FORTRAN(nident,(&ics[lprev],&inode1,&ilength,&id));
291  if(id>0){
292  if(ics[lprev+id-1]==inode1){
293  icomplex1=ij+1;
294  break;
295  }
296  }
297  }
298  }
299 // idof1=nactdof[mt*inode1+nodempc[3*index1-2]-4];
300  idof1=nactdof[mt*(inode1-1)+nodempc[3*index1-2]];
301  index2=index1;
302  while(1){
303  inode2=nodempc[3*index2-3];
304  icomplex2=0;
305  if(strcmp1(&labmpc[(id1-1)*20],"CYCLIC")==0){
306  icomplex2=atoi(&labmpc[20*(id1-1)+6]);
307  }
308  else if(strcmp1(&labmpc[(id1-1)*20],"SUBCYCLIC")==0){
309  for(ij=0;ij<*mcs;ij++){
310  ilength=cs[17*ij+3];
311  lprev=cs[17*ij+13];
312  FORTRAN(nident,(&ics[lprev],&inode2,&ilength,&id));
313  if(id>0){
314  if(ics[lprev+id-1]==inode2){
315  icomplex2=ij+1;
316  break;
317  }
318  }
319  }
320  }
321 // idof2=nactdof[mt*inode2+nodempc[3*index2-2]-4];
322  idof2=nactdof[mt*(inode2-1)+nodempc[3*index2-2]];
323  if((idof1!=0)&&(idof2!=0)){
324  insert(ipointer,&mast1,&irow,&idof1,&idof2,&ifree,&nzs_);
325  kdof1=idof1+neq[0];kdof2=idof2+neq[0];
326  insert(ipointer,&mast1,&irow,&kdof1,&kdof2,&ifree,&nzs_);
327  if(((icomplex1!=0)||(icomplex2!=0))&&
328  (icomplex1!=icomplex2)){
329  /* if(((icomplex1!=0)||(icomplex2!=0))&&
330  ((icomplex1==0)||(icomplex2==0))){*/
331  insert(ipointer,&mast1,&irow,&kdof1,&idof2,&ifree,&nzs_);
332  insert(ipointer,&mast1,&irow,&idof1,&kdof2,&ifree,&nzs_);
333  }
334  }
335  index2=nodempc[3*index2-1];
336  if(index2==0) break;
337  }
338  index1=nodempc[3*index1-1];
339  if(index1==0) break;
340  }
341  }
342 
343  else{
344 
345  /* MPC id1 /MPC id2 */
346 
347  ist1=ipompc[id1-1];
348  index1=nodempc[3*ist1-1];
349  if(index1==0) continue;
350  while(1){
351  inode1=nodempc[3*index1-3];
352  icomplex1=0;
353  if(strcmp1(&labmpc[(id1-1)*20],"CYCLIC")==0){
354  icomplex1=atoi(&labmpc[20*(id1-1)+6]);
355  }
356  else if(strcmp1(&labmpc[(id1-1)*20],"SUBCYCLIC")==0){
357  for(ij=0;ij<*mcs;ij++){
358  ilength=cs[17*ij+3];
359  lprev=cs[17*ij+13];
360  FORTRAN(nident,(&ics[lprev],&inode1,&ilength,&id));
361  if(id>0){
362  if(ics[lprev+id-1]==inode1){
363  icomplex1=ij+1;
364  break;
365  }
366  }
367  }
368  }
369 // idof1=nactdof[mt*inode1+nodempc[3*index1-2]-4];
370  idof1=nactdof[mt*(inode1-1)+nodempc[3*index1-2]];
371  ist2=ipompc[id2-1];
372  index2=nodempc[3*ist2-1];
373  if(index2==0){
374  index1=nodempc[3*index1-1];
375  if(index1==0){break;}
376  else{continue;}
377  }
378  while(1){
379  inode2=nodempc[3*index2-3];
380  icomplex2=0;
381  if(strcmp1(&labmpc[(id2-1)*20],"CYCLIC")==0){
382  icomplex2=atoi(&labmpc[20*(id2-1)+6]);
383  }
384  else if(strcmp1(&labmpc[(id2-1)*20],"SUBCYCLIC")==0){
385  for(ij=0;ij<*mcs;ij++){
386  ilength=cs[17*ij+3];
387  lprev=cs[17*ij+13];
388  FORTRAN(nident,(&ics[lprev],&inode2,&ilength,&id));
389  if(id>0){
390  if(ics[lprev+id-1]==inode2){
391  icomplex2=ij+1;
392  break;
393  }
394  }
395  }
396  }
397 // idof2=nactdof[mt*inode2+nodempc[3*index2-2]-4];
398  idof2=nactdof[mt*(inode2-1)+nodempc[3*index2-2]];
399  if((idof1!=0)&&(idof2!=0)){
400  insert(ipointer,&mast1,&irow,&idof1,&idof2,&ifree,&nzs_);
401  kdof1=idof1+neq[0];kdof2=idof2+neq[0];
402  insert(ipointer,&mast1,&irow,&kdof1,&kdof2,&ifree,&nzs_);
403  if(((icomplex1!=0)||(icomplex2!=0))&&
404  (icomplex1!=icomplex2)){
405  /* if(((icomplex1!=0)||(icomplex2!=0))&&
406  ((icomplex1==0)||(icomplex2==0))){*/
407  insert(ipointer,&mast1,&irow,&kdof1,&idof2,&ifree,&nzs_);
408  insert(ipointer,&mast1,&irow,&idof1,&kdof2,&ifree,&nzs_);
409  }
410  }
411  index2=nodempc[3*index2-1];
412  if(index2==0) break;
413  }
414  index1=nodempc[3*index1-1];
415  if(index1==0) break;
416  }
417  }
418  }
419  }
420  }
421  }
422  }
423 
424  neq[0]=2*neq[0];
425  neq[1]=neq[0];
426 
427  /* ordering the nonzero nodes in the SUPERdiagonal columns
428  mast1 contains the row numbers column per column,
429  irow the column numbers */
430 
431 /* for(i=0;i<neq[0];++i){
432  itot=0;
433  if(ipointer[i]==0){
434  printf("*ERROR in mastructcs: zero column");
435  FORTRAN(stop,());
436  }
437  istart=ipointer[i];
438  while(1){
439  ++itot;
440  ikcol[itot-1]=mast1[istart-1];
441  istart=irow[istart-1];
442  if(istart==0) break;
443  }
444  FORTRAN(isortii,(ikcol,icol,&itot,&kflag));
445  istart=ipointer[i];
446  for(j=0;j<itot-1;++j){
447  mast1[istart-1]=ikcol[j];
448  istartold=istart;
449  istart=irow[istart-1];
450  irow[istartold-1]=i+1;
451  }
452  mast1[istart-1]=ikcol[itot-1];
453  irow[istart-1]=i+1;
454  }*/
455 
456 
457  for(i=0;i<neq[0];++i){
458  if(ipointer[i]==0){
459  if(i>=neq[1]) continue;
460  printf("*ERROR in mastructcs: zero column\n");
461  FORTRAN(stop,());
462  }
463  istart=ipointer[i];
464  while(1){
465  istartold=istart;
466  istart=irow[istart-1];
467  irow[istartold-1]=i+1;
468  if(istart==0) break;
469  }
470  }
471 
472  if(neq[0]==0){
473  printf("\n*WARNING: no degrees of freedom in the model\n");
474  FORTRAN(stop,());
475  }
476 
477  printf(" number of equations\n");
478  printf(" %" ITGFORMAT "\n",neq[0]);
479  printf(" number of nonzero lower triangular matrix elements\n");
480  printf(" %" ITGFORMAT "\n",ifree-neq[0]);
481 
482  /* new meaning of icol,j1,mast1,irow:
483 
484  - irow is going to contain the row numbers of the SUBdiagonal
485  nonzero's, column per column
486  - mast1 contains the column numbers
487  - icol(i)=# SUBdiagonal nonzero's in column i
488  - jq(i)= location in field irow of the first SUBdiagonal
489  nonzero in column i
490 
491  */
492 
493  nmast=ifree;
494 
495  /* switching from a SUPERdiagonal inventory to a SUBdiagonal one */
496 
497  FORTRAN(isortii,(mast1,irow,&nmast,&kflag));
498 
499  /* filtering out the diagonal elements and generating icol and jq */
500 
501  isubtract=0;
502  for(i=0;i<neq[0];++i){icol[i]=0;}
503  k=0;
504  for(i=0;i<nmast;++i){
505  if(mast1[i]==irow[i]){++isubtract;}
506  else{
507  mast1[i-isubtract]=mast1[i];
508  irow[i-isubtract]=irow[i];
509  if(k!=mast1[i]){
510  for(l=k;l<mast1[i];++l){jq[l]=i+1-isubtract;}
511  k=mast1[i];
512  }
513  ++icol[k-1];
514  }
515  }
516  nmast=nmast-isubtract;
517  for(l=k;l<neq[0]+1;++l){jq[l]=nmast+1;}
518 
519  for(i=0;i<neq[0];++i){
520  if(jq[i+1]-jq[i]>0){
521  isize=jq[i+1]-jq[i];
522  FORTRAN(isortii,(&irow[jq[i]-1],&mast1[jq[i]-1],&isize,&kflag));
523  }
524  }
525 
526  nzs[0]=jq[neq[0]-1]-1;
527  nzs[1]=nzs[0];
528  nzs[2]=nzs[0];
529 
530  *mast1p=mast1;
531  *irowp=irow;
532 
533  return;
534 
535 }
#define ITGFORMAT
Definition: CalculiX.h:52
void FORTRAN(addimdnodecload,(ITG *nodeforc, ITG *i, ITG *imdnode, ITG *nmdnode, double *xforc, ITG *ikmpc, ITG *ilmpc, ITG *ipompc, ITG *nodempc, ITG *nmpc, ITG *imddof, ITG *nmddof, ITG *nactdof, ITG *mi, ITG *imdmpc, ITG *nmdmpc, ITG *imdboun, ITG *nmdboun, ITG *ikboun, ITG *nboun, ITG *ilboun, ITG *ithermal))
ITG strcmp1(const char *s1, const char *s2)
Definition: strcmp1.c:24
void insert(ITG *ipointer, ITG **mast1p, ITG **mast2p, ITG *i1, ITG *i2, ITG *ifree, ITG *nzs_)
Definition: insert.c:24
subroutine stop()
Definition: stop.f:19
subroutine isortii(ix, iy, n, kflag)
Definition: isortii.f:5
subroutine nident(x, px, n, id)
Definition: nident.f:25
#define ITG
Definition: CalculiX.h:51
Hosted by OpenAircraft.com, (Michigan UAV, LLC)