1 /* CalculiX - A 3-dimensional finite element program */
2 /* Copyright (C) 1998-2021 Guido Dhondt */
3
4 /* This program is free software; you can redistribute it and/or */
5 /* modify it under the terms of the GNU General Public License as */
6 /* published by the Free Software Foundation(version 2); */
7 /* */
8
9 /* This program is distributed in the hope that it will be useful, */
10 /* but WITHOUT ANY WARRANTY; without even the implied warranty of */
11 /* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the */
12 /* GNU General Public License for more details. */
13
14 /* You should have received a copy of the GNU General Public License */
15 /* along with this program; if not, write to the Free Software */
16 /* Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. */
17
18 #include <stdio.h>
19 #include <math.h>
20 #include <stdlib.h>
21 #include <string.h>
22 #include "CalculiX.h"
23
radcyc(ITG * nk,ITG * kon,ITG * ipkon,char * lakon,ITG * ne,double * cs,ITG * mcs,ITG * nkon,ITG * ialset,ITG * istartset,ITG * iendset,ITG ** kontrip,ITG * ntri,double ** cop,double ** voldp,ITG * ntrit,ITG * inocs,ITG * mi)24 void radcyc(ITG *nk,ITG *kon,ITG *ipkon,char *lakon,ITG *ne,
25 double *cs, ITG *mcs, ITG *nkon,ITG *ialset, ITG *istartset,
26 ITG *iendset,ITG **kontrip,ITG *ntri,
27 double **cop, double **voldp,ITG *ntrit, ITG *inocs,
28 ITG *mi){
29
30 /* duplicates triangular faces for cyclic radiation conditions */
31
32 char *filab=NULL;
33
34 ITG i,is,nsegments,idtie,nkt,icntrl,imag=0,*kontri=NULL,mt=mi[1]+1,
35 node,i1,i2,nope,iel,indexe,j,k,ielset,node1,node2,node3,l,jj;
36
37 double *vt=NULL,*fnt=NULL,*stnt=NULL,*eent=NULL,*qfnt=NULL,t[3],theta,
38 pi,*v=NULL,*fn=NULL,*stn=NULL,*een=NULL,*qfn=NULL,*co=NULL,
39 *vold=NULL,*emnt=NULL,*emn=NULL;
40
41 pi=4.*atan(1.);
42
43 kontri=*kontrip;co=*cop;vold=*voldp;
44
45 /* determining the maximum number of sectors */
46
47 nsegments=1;
48 for(j=0;j<*mcs;j++){
49 if(cs[17*j]>nsegments) nsegments=(ITG)(cs[17*j]);
50 }
51
52 /* assigning nodes and elements to sectors */
53
54 ielset=cs[12];
55 if((*mcs!=1)||(ielset!=0)){
56 for(i=0;i<*nk;i++) inocs[i]=-1;
57 }
58
59 for(i=0;i<*mcs;i++){
60 is=cs[17*i+4];
61 if(is==1) continue;
62 ielset=cs[17*i+12];
63 if(ielset==0) continue;
64 for(i1=istartset[ielset-1]-1;i1<iendset[ielset-1];i1++){
65 if(ialset[i1]>0){
66 iel=ialset[i1]-1;
67 if(ipkon[iel]<0) continue;
68 indexe=ipkon[iel];
69 if(strcmp1(&lakon[8*iel+3],"2")==0)nope=20;
70 else if (strcmp1(&lakon[8*iel+3],"8")==0)nope=8;
71 else if (strcmp1(&lakon[8*iel+3],"10")==0)nope=10;
72 else if (strcmp1(&lakon[8*iel+3],"4")==0)nope=4;
73 else if (strcmp1(&lakon[8*iel+3],"15")==0)nope=15;
74 else {nope=6;}
75 for(i2=0;i2<nope;++i2){
76 node=kon[indexe+i2]-1;
77 inocs[node]=i;
78 }
79 }
80 else{
81 iel=ialset[i1-2]-1;
82 do{
83 iel=iel-ialset[i1];
84 if(iel>=ialset[i1-1]-1) break;
85 if(ipkon[iel]<0) continue;
86 indexe=ipkon[iel];
87 if(strcmp1(&lakon[8*iel+3],"2")==0)nope=20;
88 else if (strcmp1(&lakon[8*iel+3],"8")==0)nope=8;
89 else if (strcmp1(&lakon[8*iel+3],"10")==0)nope=10;
90 else if (strcmp1(&lakon[8*iel+3],"4")==0)nope=4;
91 else if (strcmp1(&lakon[8*iel+3],"15")==0)nope=15;
92 else {nope=6;}
93 for(i2=0;i2<nope;++i2){
94 node=kon[indexe+i2]-1;
95 inocs[node]=i;
96 }
97 }while(1);
98 }
99 }
100 }
101
102 /* duplicating triangular faces
103 only those faces are duplicated the nodes of which belong to
104 the same cyclic symmetry. non-integer cyclic symmety numbers are
105 reduced to the next lower integer. */
106
107 *ntrit=nsegments**ntri;
108 RENEW(kontri,ITG,4**ntrit);
109 for(i=4**ntri;i<4**ntrit;i++) kontri[i]=0;
110
111 for(i=0;i<*ntri;i++){
112 node1=kontri[4*i];
113 if(inocs[node1-1]<0) continue;
114 idtie=inocs[node1-1];
115 node2=kontri[4*i+1];
116 if((inocs[node2-1]<0)||(inocs[node2-1]!=idtie)) continue;
117 node3=kontri[4*i+2];
118 if((inocs[node3-1]<0)||(inocs[node3-1]!=idtie)) continue;
119 idtie=cs[17*idtie];
120 for(k=1;k<idtie;k++){
121 j=i+k**ntri;
122 kontri[4*j]=node1+k**nk;
123 kontri[4*j+1]=node2+k**nk;
124 kontri[4*j+2]=node3+k**nk;
125 kontri[4*j+3]=kontri[4*i+3];
126 }
127 }
128
129 RENEW(co,double,3**nk*nsegments);
130 RENEW(vold,double,mt**nk*nsegments);
131 nkt=*nk*nsegments;
132
133 /* generating the coordinates for the other sectors */
134
135 icntrl=1;
136
137 FORTRAN(rectcyl,(co,v,fn,stn,qfn,een,cs,nk,&icntrl,t,filab,&imag,mi,emn));
138
139 for(jj=0;jj<*mcs;jj++){
140 is=(ITG)(cs[17*jj]);
141 for(i=1;i<is;i++){
142
143 theta=i*2.*pi/cs[17*jj];
144
145 for(l=0;l<*nk;l++){
146 if(inocs[l]==jj){
147 co[3*l+i*3**nk]=co[3*l];
148 co[1+3*l+i*3**nk]=co[1+3*l]-theta;
149 co[2+3*l+i*3**nk]=co[2+3*l];
150 }
151 }
152 }
153 }
154
155 icntrl=-1;
156
157 FORTRAN(rectcyl,(co,vt,fnt,stnt,qfnt,eent,cs,&nkt,&icntrl,t,filab,
158 &imag,mi,emnt));
159
160 *kontrip=kontri;*cop=co;*voldp=vold;
161
162 return;
163 }
164
165