CalculiX  2.13
A Free Software Three-Dimensional Structural Finite Element Program
gen3delem.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine gen3delem (kon, ipkon, lakon, ne, ipompc, nodempc, coefmpc, nmpc, nmpc_, mpcfree, ikmpc, ilmpc, labmpc, ikboun, ilboun, nboun, nboun_, nodeboun, ndirboun, xboun, iamboun, nam, inotr, trab, nk, nk_, iponoel, inoel, iponor, xnor, thicke, thickn, knor, istep, offset, t0, t1, ikforc, ilforc, rig, nforc, nforc_, nodeforc, ndirforc, xforc, iamforc, sideload, nload, ithermal, ntrans, co, ixfree, ikfree, inoelfree, iponoelmax, iperturb, tinc, tper, tmin, tmax, ctrl, typeboun, nmethod, nset, set, istartset, iendset, ialset, prop, ielprop, vold, mi, nkon, ielmat, icomposite, t0g, t1g, idefforc, iamt1, orname, orab, norien, norien_, ielorien)
 

Function/Subroutine Documentation

◆ gen3delem()

subroutine gen3delem ( integer, dimension(*)  kon,
integer, dimension(*)  ipkon,
character*8, dimension(*)  lakon,
integer  ne,
integer, dimension(*)  ipompc,
integer, dimension(3,*)  nodempc,
real*8, dimension(*)  coefmpc,
integer  nmpc,
integer  nmpc_,
integer  mpcfree,
integer, dimension(*)  ikmpc,
integer, dimension(*)  ilmpc,
character*20, dimension(*)  labmpc,
integer, dimension(*)  ikboun,
integer, dimension(*)  ilboun,
integer  nboun,
integer  nboun_,
integer, dimension(*)  nodeboun,
integer, dimension(*)  ndirboun,
real*8, dimension(*)  xboun,
integer, dimension(*)  iamboun,
integer  nam,
integer, dimension(2,*)  inotr,
real*8, dimension(7,*)  trab,
integer  nk,
integer  nk_,
integer, dimension(*)  iponoel,
integer, dimension(3,*)  inoel,
integer, dimension(2,*)  iponor,
real*8, dimension(*)  xnor,
real*8, dimension(mi(3),*)  thicke,
real*8, dimension(2,*)  thickn,
integer, dimension(*)  knor,
integer  istep,
real*8, dimension(2,*)  offset,
real*8, dimension(*)  t0,
real*8, dimension(*)  t1,
integer, dimension(*)  ikforc,
integer, dimension(*)  ilforc,
integer, dimension(*)  rig,
integer  nforc,
integer  nforc_,
integer, dimension(2,*)  nodeforc,
integer, dimension(*)  ndirforc,
real*8, dimension(*)  xforc,
integer, dimension(*)  iamforc,
character*20, dimension(*)  sideload,
integer  nload,
integer, dimension(2)  ithermal,
integer  ntrans,
real*8, dimension(3,*)  co,
integer  ixfree,
integer  ikfree,
integer  inoelfree,
integer  iponoelmax,
integer  iperturb,
real*8  tinc,
real*8  tper,
real*8  tmin,
real*8  tmax,
real*8, dimension(*)  ctrl,
character*1, dimension(*)  typeboun,
integer  nmethod,
integer  nset,
character*81, dimension(*)  set,
integer, dimension(*)  istartset,
integer, dimension(*)  iendset,
integer, dimension(*)  ialset,
real*8, dimension(*)  prop,
integer, dimension(*)  ielprop,
real*8, dimension(0:mi(2),*)  vold,
integer, dimension(*)  mi,
integer  nkon,
integer, dimension(mi(3),*)  ielmat,
integer  icomposite,
real*8, dimension(2,*)  t0g,
real*8, dimension(2,*)  t1g,
integer, dimension(*)  idefforc,
integer, dimension(*)  iamt1,
character*80, dimension(*)  orname,
real*8, dimension(7,*)  orab,
integer  norien,
integer  norien_,
integer, dimension(mi(3),*)  ielorien 
)
30 !
31 ! generates three-dimensional elements:
32 ! for plane stress
33 ! for plane strain
34 ! for plate and shell elements
35 ! for beam elements
36 !
37  implicit none
38 !
39  character*1 typeboun(*)
40  character*8 lakon(*)
41  character*20 labmpc(*),sideload(*)
42  character*80 orname(*)
43  character*81 set(*)
44 !
45  integer ipompc(*),nodempc(3,*),nmpc,nmpc_,mpcfree,ikmpc(*),
46  & ilmpc(*),kon(*),ipkon(*),ne,indexe,i,j,node,
47  & ikboun(*),ilboun(*),nboun,nboun_,
48  & neigh(7,8),nodeboun(*),ndirboun(*),nk,
49  & nk_,iponoel(*),inoel(3,*),inoelfree,istep,nmpcold,
50  & ikforc(*),ilforc(*),nodeforc(2,*),ndirforc(*),iamforc(*),
51  & nforc,nforc_,ithermal(2),nload,iamboun(*),
52  & ntrans,inotr(2,*),nam,iponoelmax,iperturb,numnod,itransaxial,
53  & rig(*),nmethod,nset,istartset(*),iendset(*),ialset(*),nkon,
54  & ielprop(*),mi(*),nope,
55  & ielmat(mi(3),*),iponor(2,*),knor(*),ixfree,ikfree,icomposite,
56  & idefforc(*),idim,iamt1(*),norien,norien_,ielorien(mi(3),*)
57 !
58  real*8 coefmpc(*),thicke(mi(3),*),xnor(*),thickn(2,*),tinc,
59  & tper,tmin,t0g(2,*),t1g(2,*),e1(3),e2(3),xt1(3),orab(7,*),
60  & tmax,offset(2,*),t0(*),t1(*),xforc(*),trab(7,*),co(3,*),
61  & xboun(*),pi,ctrl(*),prop(*),vold(0:mi(2),*),
62  & coloc(3,8)
63 !
64  data neigh /1,9,2,12,4,17,5,2,9,1,10,3,18,6,
65  & 3,11,4,10,2,19,7,4,11,3,12,1,20,8,
66  & 5,13,6,16,8,17,1,6,13,5,14,7,18,2,
67  & 7,15,8,14,6,19,3,8,15,7,16,5,20,4/
68 !
69  data coloc /-1.,-1.,-1.,1.,-1.,-1.,1.,1.,-1.,-1.,1.,-1.,
70  & -1.,-1.,1.,1.,-1.,1.,1.,1.,1.,-1.,1.,1./
71 !
72  pi=4.d0*datan(1.d0)
73 !
74 ! catalogueing the element per node relationship for shell/beam
75 ! elements and transferring the nodal thickness to the elements
76 !
77 ! inoelfree=1 means that there is at least one 1D or 2D element
78 ! in the structure. Otherwise inoelfree=0.
79 !
80  if((istep.eq.1).and.(inoelfree.eq.1)) then
81 !
82 ! shift of the connectivity for composite elements
83 !
84  if(icomposite.eq.1) then
85  call changekon(ne,ipkon,lakon,mi,nkon,thicke,ielmat,kon)
86  endif
87 !
88  itransaxial=0
89 !
90  do i=1,ne
91  if(ipkon(i).lt.0) cycle
92  if((lakon(i)(1:2).ne.'C3').and.(lakon(i)(1:1).ne.'D').and.
93  & (lakon(i)(1:1).ne.'G').and.(lakon(i)(1:1).ne.'E').and.
94  & (lakon(i)(1:4).ne.'MASS')) then
95 !
96 ! number of nodes belonging to the element
97 !
98  if((lakon(i)(1:3).eq.'B31').or.
99  & (lakon(i)(1:4).eq.'T3D2')) then
100  numnod=2
101  elseif((lakon(i)(1:1).eq.'B').or.
102  & (lakon(i)(1:1).eq.'T')) then
103  numnod=3
104  elseif((lakon(i)(2:2).eq.'3').or.
105  & (lakon(i)(4:4).eq.'3')) then
106  numnod=3
107  elseif((lakon(i)(2:2).eq.'4').or.
108  & (lakon(i)(4:4).eq.'4')) then
109  numnod=4
110  elseif((lakon(i)(2:2).eq.'6').or.
111  & (lakon(i)(4:4).eq.'6')) then
112  numnod=6
113  elseif((lakon(i)(2:2).eq.'8').or.
114  & (lakon(i)(4:4).eq.'8')) then
115  numnod=8
116  endif
117 !
118  indexe=ipkon(i)
119  do j=1,numnod
120  node=kon(indexe+j)
121  iponoelmax=max(iponoelmax,node)
122  inoel(1,inoelfree)=i
123  inoel(2,inoelfree)=j
124  inoel(3,inoelfree)=iponoel(node)
125  iponoel(node)=inoelfree
126  inoelfree=inoelfree+1
127  if(lakon(i)(1:2).ne.'CA') then
128  if(thickn(1,node).gt.0.d0)
129  & thicke(1,indexe+j)=thickn(1,node)
130  if(mi(3).gt.1) then
131  if(thickn(2,node).gt.0.d0)
132  & thicke(2,indexe+j)=thickn(2,node)
133  endif
134  endif
135  if(thicke(1,indexe+j).le.0.d0) then
136  if(lakon(i)(1:1).eq.'C') then
137 !
138 ! default for plane stress and plane strain elements
139 !
140  thicke(1,indexe+j)=1.d0
141  else
142  write(*,*)'*ERROR in gen3delem: first thickness'
143  write(*,*)' in node ',j,' of element ',i
144  write(*,*)' is zero'
145  call exit(201)
146  endif
147  endif
148  if((lakon(i)(1:1).eq.'B').and.
149  & (thicke(2,indexe+j).le.0.d0)) then
150  write(*,*) '*ERROR in gen3delem: second thickness'
151  write(*,*)' in node ',j,' of beam element ',i
152  write(*,*)' is zero'
153  call exit(201)
154  endif
155  enddo
156  endif
157  enddo
158 !
159 ! calculating the normals in nodes belonging to shells/beams
160 !
161  nmpcold=nmpc
162 !
163  call gen3dnor(nk,nk_,co,iponoel,inoel,iponoelmax,kon,ipkon,
164  & lakon,ne,thicke,offset,iponor,xnor,knor,rig,iperturb,tinc,
165  & tper,tmin,tmax,ctrl,ipompc,nodempc,coefmpc,nmpc,nmpc_,
166  & mpcfree,ikmpc,ilmpc,labmpc,ikboun,ilboun,nboun,nboun_,
167  & nodeboun,ndirboun,xboun,iamboun,typeboun,nam,ntrans,inotr,
168  & trab,ikfree,ixfree,nmethod,ithermal,istep,mi,icomposite,
169  & ielmat,vold)
170 !
171  endif
172 !
173  if(istep.eq.1) then
174 !
175 ! if there is any plane stress, plane strain or axisymmetric
176 ! element the structure should lie in the z=0 plane
177 !
178  if(inoelfree.ne.0) then
179  do i=1,ne
180  if(ipkon(i).lt.0) cycle
181  if((lakon(i)(1:2).eq.'CP').or.
182  & (lakon(i)(1:2).eq.'CA')) then
183  indexe=ipkon(i)
184  read(lakon(i)(4:4),'(i1)') nope
185  do j=1,nope
186  node=kon(indexe+j)
187  if(dabs(co(3,node)).gt.0.d0) then
188  write(*,*) '*ERROR in gen3delem. The structure'
189  write(*,*) ' contains plane stress, plane'
190  write(*,*) ' strain or axisymmetric'
191  write(*,*) ' elements and should lie in '
192  write(*,*) ' the z=0 plane. This is at'
193  write(*,*) ' least not the case for node',
194  & node
195  call exit(201)
196  endif
197  enddo
198  endif
199  enddo
200  endif
201 !
202 ! 1D and 2D elements
203 !
204  if(inoelfree.ne.0) then
205  do i=1,ne
206  if(ipkon(i).lt.0) cycle
207  if((lakon(i)(1:2).eq.'CP').or.
208  & (lakon(i)(1:1).eq.'S').or.
209  & (lakon(i)(1:2).eq.'M3').or.
210  & (lakon(i)(1:2).eq.'CA')) then
211 !
212  call gen3dfrom2d(i,kon,ipkon,lakon,ne,iponor,xnor,knor,
213  & thicke,offset,ntrans,inotr,trab,ikboun,ilboun,nboun,
214  & nboun_,nodeboun,ndirboun,xboun,iamboun,typeboun,ipompc,
215  & nodempc,coefmpc,nmpc,nmpc_,mpcfree,ikmpc,ilmpc,labmpc,
216  & nk,nk_,co,rig,nmethod,iperturb,ithermal,mi,nam,
217  & icomposite,ielmat,vold,orname,orab,norien,norien_,
218  & ielorien)
219 !
220  elseif((lakon(i)(1:1).eq.'B').or.
221  & (lakon(i)(1:1).eq.'T')) then
222  call gen3dfrom1d(i,kon,ipkon,lakon,ne,iponor,xnor,knor,
223  & thicke,ntrans,inotr,trab,nk,nk_,co,offset,mi)
224  endif
225 !
226  if(lakon(i)(1:4).eq.'CPE3') then
227  lakon(i)(1:7)='C3D6 E'
228  elseif(lakon(i)(1:5).eq.'CPE4R') then
229  lakon(i)(1:7)='C3D8R E'
230  elseif(lakon(i)(1:4).eq.'CPE4') then
231  lakon(i)(1:7)='C3D8 E'
232  elseif(lakon(i)(1:4).eq.'CPE6') then
233  lakon(i)(1:7)='C3D15 E'
234  elseif(lakon(i)(1:5).eq.'CPE8R') then
235  lakon(i)(1:7)='C3D20RE'
236  elseif(lakon(i)(1:4).eq.'CPE8') then
237  lakon(i)(1:7)='C3D20 E'
238  elseif(lakon(i)(1:4).eq.'CPS3') then
239  lakon(i)(1:7)='C3D6 S'
240  elseif(lakon(i)(1:5).eq.'CPS4R') then
241  lakon(i)(1:7)='C3D8R S'
242  elseif(lakon(i)(1:4).eq.'CPS4') then
243  lakon(i)(1:7)='C3D8 S'
244  elseif(lakon(i)(1:4).eq.'CPS6') then
245  lakon(i)(1:7)='C3D15 S'
246  elseif(lakon(i)(1:5).eq.'CPS8R') then
247  lakon(i)(1:7)='C3D20RS'
248  elseif(lakon(i)(1:4).eq.'CPS8') then
249  lakon(i)(1:7)='C3D20 S'
250  elseif(lakon(i)(1:4).eq.'CAX3') then
251  lakon(i)(1:7)='C3D6 A'
252  elseif(lakon(i)(1:5).eq.'CAX4R') then
253  lakon(i)(1:7)='C3D8R A'
254  elseif(lakon(i)(1:4).eq.'CAX4') then
255  lakon(i)(1:7)='C3D8 A'
256  elseif(lakon(i)(1:4).eq.'CAX6') then
257  lakon(i)(1:7)='C3D15 A'
258  elseif(lakon(i)(1:5).eq.'CAX8R') then
259  lakon(i)(1:7)='C3D20RA'
260  elseif(lakon(i)(1:4).eq.'CAX8') then
261  lakon(i)(1:7)='C3D20 A'
262  elseif((lakon(i)(1:2).eq.'S3').or.
263  & (lakon(i)(1:4).eq.'M3D3')) then
264  lakon(i)(1:7)='C3D6 L'
265  elseif((lakon(i)(1:3).eq.'S4R').or.
266  & (lakon(i)(1:5).eq.'M3D4R')) then
267  lakon(i)(1:7)='C3D8R L'
268  elseif((lakon(i)(1:2).eq.'S4').or.
269  & (lakon(i)(1:4).eq.'M3D4')) then
270  lakon(i)(1:7)='C3D8I L'
271  elseif((lakon(i)(1:2).eq.'S6').or.
272  & (lakon(i)(1:4).eq.'M3D6')) then
273  lakon(i)(1:7)='C3D15 L'
274  elseif((lakon(i)(1:3).eq.'S8R').or.
275  & (lakon(i)(1:5).eq.'M3D8R')) then
276  lakon(i)(1:7)='C3D20RL'
277  elseif((lakon(i)(1:2).eq.'S8').or.
278  & (lakon(i)(1:4).eq.'M3D8')) then
279  lakon(i)(1:7)='C3D20 L'
280  elseif(lakon(i)(1:4).eq.'B31R') then
281  lakon(i)(1:7)='C3D8R B'
282  elseif((lakon(i)(1:3).eq.'B31').or.
283  & (lakon(i)(1:4).eq.'T3D2')) then
284  lakon(i)(1:7)='C3D8I B'
285  elseif(lakon(i)(1:4).eq.'B32R') then
286  lakon(i)(1:7)='C3D20RB'
287  elseif((lakon(i)(1:3).eq.'B32').or.
288  & (lakon(i)(1:4).eq.'T3D3')) then
289  lakon(i)(1:7)='C3D20 B'
290  endif
291  enddo
292 c Bernhardi start
293  endif
294  do i=1,ne
295  if(lakon(i)(1:5).eq.'C3D8I') then
296  call genmodes(i,kon,ipkon,lakon,ne,nk,nk_,co)
297  endif
298  enddo
299 c Bernhardi end
300 !
301 ! filling the new KNOT MPC's (needs the coordinates
302 ! of the expanded nodes)
303 !
304  if(inoelfree.ne.0) then
305  call fillknotmpc(co,ipompc,nodempc,coefmpc,labmpc,
306  & nmpc,nmpcold,mpcfree,idim,e1,e2,xt1)
307  endif
308 !
309  endif
310 !
311 ! generating MPC's to connect shells and beams with solid
312 ! elements or spring elements or mass elements
313 !
314  if((inoelfree.ne.0).and.(istep.eq.1)) then
315  call gen3dconnect(kon,ipkon,lakon,ne,iponoel,inoel,
316  & iponoelmax,rig,iponor,xnor,knor,ipompc,nodempc,coefmpc,nmpc,
317  & nmpc_,mpcfree,ikmpc,ilmpc,labmpc)
318  endif
319 !
320  if(inoelfree.ne.0) then
321 !
322 ! multiplying existing boundary conditions
323 !
324  call gen3dboun(ikboun,ilboun,nboun,nboun_,nodeboun,ndirboun,
325  & xboun,iamboun,typeboun,iponoel,inoel,iponoelmax,kon,ipkon,
326  & lakon,ne,iponor,xnor,knor,ipompc,nodempc,coefmpc,nmpc,nmpc_,
327  & mpcfree,ikmpc,ilmpc,labmpc,rig,ntrans,inotr,trab,nam,nk,nk_,
328  & co,nmethod,iperturb,istep,vold,mi)
329 !
330 ! updating the nodal surfaces: establishing links between the user
331 ! defined nodes and the newly generated nodes (mid-nodes
332 ! for 2d elements, mean of corner nodes for 1d elements)
333 !
334  if(istep.eq.1) then
335  call gen3dsurf(iponoel,inoel,iponoelmax,kon,ipkon,
336  & lakon,ne,iponor,knor,ipompc,nodempc,coefmpc,nmpc,nmpc_,
337  & mpcfree,ikmpc,ilmpc,labmpc,rig,ntrans,inotr,trab,nam,nk,
338  & nk_,co,nmethod,iperturb,nset,set,istartset,iendset,ialset,
339  & ikboun,ilboun,nboun,nboun_,nodeboun,ndirboun,xboun,
340  & iamboun,typeboun,mi,vold)
341  endif
342 !
343 ! updating the MPCs: establishing links between the user
344 ! defined nodes and the newly generated nodes (mid-nodes
345 ! for 2d elements, mean of corner nodes for 1d elements)
346 !
347 ! is needed in each step since new SPC's in local
348 ! coordinates can be defined which correspond to
349 ! MPC's in global coordinates
350 !
351 c if(istep.eq.1) then
352  call gen3dmpc(ipompc,nodempc,coefmpc,nmpc,nmpc_,mpcfree,
353  & ikmpc,ilmpc,labmpc,iponoel,inoel,iponoelmax,kon,ipkon,
354  & lakon,ne,iponor,xnor,knor,rig)
355 c endif
356 !
357 ! updating the temperatures
358 !
359  if(ithermal(1).gt.0) then
360  call gen3dtemp(iponoel,inoel,iponoelmax,kon,ipkon,lakon,ne,
361  & iponor,xnor,knor,t0,t1,thicke,offset,rig,nk,nk_,co,
362  & istep,ithermal,vold,mi,t0g,t1g,nam,iamt1)
363  endif
364 !
365 ! updating the concentrated loading
366 !
367  call gen3dforc(ikforc,ilforc,nforc,nforc_,nodeforc,
368  & ndirforc,xforc,iamforc,ntrans,inotr,trab,rig,ipompc,nodempc,
369  & coefmpc,nmpc,nmpc_,mpcfree,ikmpc,ilmpc,labmpc,iponoel,inoel,
370  & iponoelmax,kon,ipkon,lakon,ne,iponor,xnor,knor,nam,nk,nk_,
371  & co,thicke,nodeboun,ndirboun,ikboun,ilboun,nboun,nboun_,
372  & iamboun,typeboun,xboun,nmethod,iperturb,istep,vold,mi,
373  & idefforc)
374  endif
375 !
376  return
subroutine gen3dforc(ikforc, ilforc, nforc, nforc_, nodeforc, ndirforc, xforc, iamforc, ntrans, inotr, trab, rig, ipompc, nodempc, coefmpc, nmpc, nmpc_, mpcfree, ikmpc, ilmpc, labmpc, iponoel, inoel, iponoelmax, kon, ipkon, lakon, ne, iponor, xnor, knor, nam, nk, nk_, co, thicke, nodeboun, ndirboun, ikboun, ilboun, nboun, nboun_, iamboun, typeboun, xboun, nmethod, iperturb, istep, vold, mi, idefforc)
Definition: gen3dforc.f:25
#define max(a, b)
Definition: cascade.c:32
subroutine gen3dconnect(kon, ipkon, lakon, ne, iponoel, inoel, iponoelmax, rig, iponor, xnor, knor, ipompc, nodempc, coefmpc, nmpc, nmpc_, mpcfree, ikmpc, ilmpc, labmpc)
Definition: gen3dconnect.f:22
subroutine gen3dsurf(iponoel, inoel, iponoelmax, kon, ipkon, lakon, ne, iponor, knor, ipompc, nodempc, coefmpc, nmpc, nmpc_, mpcfree, ikmpc, ilmpc, labmpc, rig, ntrans, inotr, trab, nam, nk, nk_, co, nmethod, iperturb, nset, set, istartset, iendset, ialset, ikboun, ilboun, nboun, nboun_, nodeboun, ndirboun, xboun, iamboun, typeboun, mi, vold)
Definition: gen3dsurf.f:25
subroutine genmodes(i, kon, ipkon, lakon, ne, nk, nk_, co)
Definition: genmodes.f:21
subroutine gen3dboun(ikboun, ilboun, nboun, nboun_, nodeboun, ndirboun, xboun, iamboun, typeboun, iponoel, inoel, iponoelmax, kon, ipkon, lakon, ne, iponor, xnor, knor, ipompc, nodempc, coefmpc, nmpc, nmpc_, mpcfree, ikmpc, ilmpc, labmpc, rig, ntrans, inotr, trab, nam, nk, nk_, co, nmethod, iperturb, istep, vold, mi)
Definition: gen3dboun.f:24
subroutine changekon(ne, ipkon, lakon, mi, nkon, thicke, ielmat, kon)
Definition: changekon.f:20
subroutine gen3dtemp(iponoel, inoel, iponoelmax, kon, ipkon, lakon, ne, iponor, xnor, knor, t0, t1, thicke, offset, rig, nk, nk_, co, istep, ithermal, vold, mi, t0g, t1g, nam, iamt1)
Definition: gen3dtemp.f:22
subroutine fillknotmpc(co, ipompc, nodempc, coefmpc, labmpc, nmpc, nmpcold, mpcfree, idim, e1, e2, t1)
Definition: fillknotmpc.f:21
subroutine gen3dnor(nk, nk_, co, iponoel, inoel, iponoelmax, kon, ipkon, lakon, ne, thicke, offset, iponor, xnor, knor, rig, iperturb, tinc, tper, tmin, tmax, ctrl, ipompc, nodempc, coefmpc, nmpc, nmpc_, mpcfree, ikmpc, ilmpc, labmpc, ikboun, ilboun, nboun, nboun_, nodeboun, ndirboun, xboun, iamboun, typeboun, nam, ntrans, inotr, trab, ikfree, ixfree, nmethod, ithermal, istep, mi, icomposite, ielmat, vold)
Definition: gen3dnor.f:25
subroutine gen3dfrom1d(i, kon, ipkon, lakon, ne, iponor, xnor, knor, thicke, ntrans, inotr, trab, nk, nk_, co, offset, mi)
Definition: gen3dfrom1d.f:21
subroutine gen3dmpc(ipompc, nodempc, coefmpc, nmpc, nmpc_, mpcfree, ikmpc, ilmpc, labmpc, iponoel, inoel, iponoelmax, kon, ipkon, lakon, ne, iponor, xnor, knor, rig)
Definition: gen3dmpc.f:22
subroutine gen3dfrom2d(i, kon, ipkon, lakon, ne, iponor, xnor, knor, thicke, offset, ntrans, inotr, trab, ikboun, ilboun, nboun, nboun_, nodeboun, ndirboun, xboun, iamboun, typeboun, ipompc, nodempc, coefmpc, nmpc, nmpc_, mpcfree, ikmpc, ilmpc, labmpc, nk, nk_, co, rig, nmethod, iperturb, ithermal, mi, nam, icomposite, ielmat, vold, orname, orab, norien, norien_, ielorien)
Definition: gen3dfrom2d.f:25
Hosted by OpenAircraft.com, (Michigan UAV, LLC)