root/src/frdcyc.c

/* [<][>][^][v][top][bottom][index][help] */

DEFINITIONS

This source file includes following definitions.
  1. frdcyc

   1 /*     CalculiX - A 3-dimensional finite element program                   */
   2 /*              Copyright (C) 1998-2015 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 
  24 void frdcyc(double *co,ITG *nk,ITG *kon,ITG *ipkon,char *lakon,ITG *ne,double *v,
  25             double *stn,ITG *inum,ITG *nmethod,ITG *kode,char *filab,
  26             double *een,double *t1,double *fn,double *time,double *epn,
  27             ITG *ielmat,char *matname, double *cs, ITG *mcs, ITG *nkon,
  28             double *enern, double *xstaten, ITG *nstate_, ITG *istep,
  29             ITG *iinc, ITG *iperturb, double *ener, ITG *mi, char *output,
  30             ITG *ithermal, double *qfn, ITG *ialset, ITG *istartset,
  31             ITG *iendset, double *trab, ITG *inotr, ITG *ntrans,
  32             double *orab, ITG *ielorien, ITG *norien, double *sti,
  33             double *veold, ITG *noddiam,char *set,ITG *nset, double *emn,
  34             double *thicke,char* jobnamec,ITG *ne0,double *cdn,ITG *mortar,ITG *nmat){
  35 
  36   /* duplicates fields for static cyclic symmetric calculations */
  37 
  38   char *lakont=NULL,description[13]="            ";
  39 
  40   ITG nkt,icntrl,*kont=NULL,*ipkont=NULL,*inumt=NULL,*ielmatt=NULL,net,i,l,
  41      imag=0,mode=-1,ngraph,*inocs=NULL,*ielcs=NULL,l1,l2,is,
  42       jj,node,i1,i2,nope,iel,indexe,j,ielset,*inotrt=NULL,mt=mi[1]+1,
  43       *ipneigh=NULL,*neigh=NULL,net0;
  44 
  45   double *vt=NULL,*fnt=NULL,*stnt=NULL,*eent=NULL,*cot=NULL,*t1t=NULL,
  46          *epnt=NULL,*enernt=NULL,*xstatent=NULL,theta,pi,t[3],*qfnt=NULL,
  47          *vr=NULL,*vi=NULL,*stnr=NULL,*stni=NULL,*vmax=NULL,*stnmax=NULL,
  48          *stit=NULL,*eenmax=NULL,*fnr=NULL,*fni=NULL,*emnt=NULL,*qfx=NULL,
  49          *cdnr=NULL,*cdni=NULL;
  50 
  51   pi=4.*atan(1.);
  52 
  53   /* determining the maximum number of sectors to be plotted */
  54 
  55   ngraph=1;
  56   for(j=0;j<*mcs;j++){
  57     if(cs[17*j+4]>ngraph) ngraph=cs[17*j+4];
  58   }
  59 
  60   /* assigning nodes and elements to sectors */
  61 
  62   NNEW(inocs,ITG,*nk);
  63   NNEW(ielcs,ITG,*ne);
  64   ielset=cs[12];
  65   if((*mcs!=1)||(ielset!=0)){
  66     for(i=0;i<*nk;i++) inocs[i]=-1;
  67     for(i=0;i<*ne;i++) ielcs[i]=-1;
  68   }
  69 
  70   for(i=0;i<*mcs;i++){
  71     is=cs[17*i+4];
  72     if(is==1) continue;
  73     ielset=cs[17*i+12];
  74     if(ielset==0) continue;
  75     for(i1=istartset[ielset-1]-1;i1<iendset[ielset-1];i1++){
  76       if(ialset[i1]>0){
  77         iel=ialset[i1]-1;
  78         if(ipkon[iel]<0) continue;
  79         ielcs[iel]=i;
  80         indexe=ipkon[iel];
  81         if(strcmp1(&lakon[8*iel+3],"2")==0)nope=20;
  82         else if (strcmp1(&lakon[8*iel+3],"8")==0)nope=8;
  83         else if (strcmp1(&lakon[8*iel+3],"10")==0)nope=10;
  84         else if (strcmp1(&lakon[8*iel+3],"4")==0)nope=4;
  85         else if (strcmp1(&lakon[8*iel+3],"15")==0)nope=15;
  86         else {nope=6;}
  87         for(i2=0;i2<nope;++i2){
  88           node=kon[indexe+i2]-1;
  89           inocs[node]=i;
  90         }
  91       }
  92       else{
  93         iel=ialset[i1-2]-1;
  94         do{
  95           iel=iel-ialset[i1];
  96           if(iel>=ialset[i1-1]-1) break;
  97           if(ipkon[iel]<0) continue;
  98           ielcs[iel]=i;
  99           indexe=ipkon[iel];
 100           if(strcmp1(&lakon[8*iel+3],"2")==0)nope=20;
 101           else if (strcmp1(&lakon[8*iel+3],"8")==0)nope=8;
 102           else if (strcmp1(&lakon[8*iel+3],"10")==0)nope=10;
 103           else if (strcmp1(&lakon[8*iel+3],"4")==0)nope=4;
 104           else if (strcmp1(&lakon[8*iel+3],"15")==0)nope=15;
 105           else {nope=6;}
 106           for(i2=0;i2<nope;++i2){
 107             node=kon[indexe+i2]-1;
 108             inocs[node]=i;
 109           }
 110         }while(1);
 111       }
 112     } 
 113   }
 114 
 115   NNEW(cot,double,3**nk*ngraph);
 116   if(*ntrans>0)NNEW(inotrt,ITG,2**nk*ngraph);
 117 
 118   if((strcmp1(&filab[0],"U ")==0)||
 119      ((strcmp1(&filab[87],"NT  ")==0)&&(*ithermal>=2)))
 120     NNEW(vt,double,mt**nk*ngraph);
 121   if((strcmp1(&filab[87],"NT  ")==0)&&(*ithermal<2))
 122     NNEW(t1t,double,*nk*ngraph);
 123   if((strcmp1(&filab[174],"S   ")==0)||(strcmp1(&filab[1044],"ZZS ")==0)||
 124      (strcmp1(&filab[1044],"ERR ")==0))
 125     NNEW(stnt,double,6**nk*ngraph);
 126   if(strcmp1(&filab[261],"E   ")==0)
 127     NNEW(eent,double,6**nk*ngraph);
 128   if((strcmp1(&filab[348],"RF  ")==0)||(strcmp1(&filab[783],"RFL ")==0))
 129     NNEW(fnt,double,mt**nk*ngraph);
 130   if(strcmp1(&filab[435],"PEEQ")==0)
 131     NNEW(epnt,double,*nk*ngraph);
 132   if(strcmp1(&filab[522],"ENER")==0)
 133     NNEW(enernt,double,*nk*ngraph);
 134   if(strcmp1(&filab[609],"SDV ")==0)
 135     NNEW(xstatent,double,*nstate_**nk*ngraph);
 136   if(strcmp1(&filab[696],"HFL ")==0)
 137     NNEW(qfnt,double,3**nk*ngraph);
 138   if((strcmp1(&filab[1044],"ZZS ")==0)||(strcmp1(&filab[1044],"ERR ")==0)||
 139      (strcmp1(&filab[2175],"CONT")==0))
 140     NNEW(stit,double,6*mi[0]**ne*ngraph);
 141   if(strcmp1(&filab[2697],"ME  ")==0)
 142     NNEW(emnt,double,6**nk*ngraph);
 143 
 144   /* the topology only needs duplication the first time it is
 145      stored in the frd file (*kode=1)
 146      the above two lines are not true: lakon is needed for
 147      contact information in frd.f */
 148 
 149 //  if(*kode==1){
 150     NNEW(kont,ITG,*nkon*ngraph);
 151     NNEW(ipkont,ITG,*ne*ngraph);
 152     NNEW(lakont,char,8**ne*ngraph);
 153     NNEW(ielmatt,ITG,mi[2]**ne*ngraph);
 154 //  }
 155   NNEW(inumt,ITG,*nk*ngraph);
 156   
 157   nkt=ngraph**nk;
 158   net0=(ngraph-1)**ne+(*ne0);
 159   net=ngraph**ne;
 160 
 161   /* copying the coordinates of the first sector */
 162   
 163   for(l=0;l<3**nk;l++){cot[l]=co[l];}
 164   if(*ntrans>0){for(l=0;l<*nk;l++){inotrt[2*l]=inotr[2*l];}}
 165 
 166   /* copying the topology of the first sector */
 167   
 168 //  if(*kode==1){
 169       for(l=0;l<*nkon;l++){kont[l]=kon[l];}
 170       for(l=0;l<*ne;l++){ipkont[l]=ipkon[l];}
 171       for(l=0;l<8**ne;l++){lakont[l]=lakon[l];}
 172       for(l=0;l<mi[2]**ne;l++){ielmatt[l]=ielmat[l];}
 173 //  }  
 174 
 175   /* generating the coordinates for the other sectors */
 176   
 177   icntrl=1;
 178   
 179   FORTRAN(rectcyl,(cot,v,fn,stn,qfn,een,cs,nk,&icntrl,t,filab,&imag,mi,emn));
 180   
 181   for(jj=0;jj<*mcs;jj++){
 182     is=cs[17*jj+4];
 183     for(i=1;i<is;i++){
 184       
 185       theta=i*2.*pi/cs[17*jj];
 186       
 187       for(l=0;l<*nk;l++){
 188         if(inocs[l]==jj){
 189           cot[3*l+i*3**nk]=cot[3*l];
 190           cot[1+3*l+i*3**nk]=cot[1+3*l]+theta;
 191           cot[2+3*l+i*3**nk]=cot[2+3*l];
 192         }
 193       }
 194       
 195       if(*ntrans>0){
 196           for(l=0;l<*nk;l++){
 197               if(inocs[l]==jj){
 198                   inotrt[2*l+i*2**nk]=inotrt[2*l];
 199               }
 200           }
 201       }
 202       
 203       //   if(*kode==1){
 204         
 205         for(l=0;l<*nkon;l++){kont[l+i**nkon]=kon[l]+i**nk;}
 206         for(l=0;l<*ne;l++){
 207           if(ielcs[l]==jj){
 208             if(ipkon[l]>=0){
 209               ipkont[l+i**ne]=ipkon[l]+i**nkon;
 210               ielmatt[mi[2]*(l+i**ne)]=ielmat[mi[2]*l];
 211               for(l1=0;l1<8;l1++){
 212                 l2=8*l+l1;
 213                 lakont[l2+i*8**ne]=lakon[l2];
 214               }
 215             }
 216             else ipkont[l+i**ne]=-1;
 217           }
 218         }
 219         //   }
 220     }
 221   }
 222 
 223   icntrl=-1;
 224     
 225   FORTRAN(rectcyl,(cot,vt,fnt,stnt,qfnt,eent,cs,&nkt,&icntrl,t,filab,
 226                    &imag,mi,emn));
 227   
 228   /* mapping the results to the other sectors */
 229   
 230   for(l=0;l<*nk;l++){inumt[l]=inum[l];}
 231   
 232   icntrl=2;
 233   
 234   FORTRAN(rectcyl,(co,v,fn,stn,qfn,een,cs,nk,&icntrl,t,filab,&imag,mi,emn));
 235   
 236   if((strcmp1(&filab[0],"U ")==0)||
 237      ((strcmp1(&filab[87],"NT  ")==0)&&(*ithermal>=2)))
 238     for(l=0;l<mt**nk;l++){vt[l]=v[l];};
 239   if((strcmp1(&filab[87],"NT  ")==0)&&(*ithermal<2))
 240     for(l=0;l<*nk;l++){t1t[l]=t1[l];};
 241   if(strcmp1(&filab[174],"S   ")==0)
 242     for(l=0;l<6**nk;l++){stnt[l]=stn[l];};
 243   if(strcmp1(&filab[261],"E   ")==0)
 244     for(l=0;l<6**nk;l++){eent[l]=een[l];};
 245   if((strcmp1(&filab[348],"RF  ")==0)||(strcmp1(&filab[783],"RFL ")==0))
 246     for(l=0;l<mt**nk;l++){fnt[l]=fn[l];};
 247   if(strcmp1(&filab[435],"PEEQ")==0)
 248     for(l=0;l<*nk;l++){epnt[l]=epn[l];};
 249   if(strcmp1(&filab[522],"ENER")==0)
 250     for(l=0;l<*nk;l++){enernt[l]=enern[l];};
 251   if(strcmp1(&filab[609],"SDV ")==0)
 252     for(l=0;l<*nstate_**nk;l++){xstatent[l]=xstaten[l];};
 253   if(strcmp1(&filab[696],"HFL ")==0)
 254     for(l=0;l<3**nk;l++){qfnt[l]=qfn[l];};
 255   if((strcmp1(&filab[1044],"ZZS ")==0)||(strcmp1(&filab[1044],"ERR ")==0)||
 256      (strcmp1(&filab[2175],"CONT")==0))
 257     for(l=0;l<6*mi[0]**ne;l++){stit[l]=sti[l];};
 258   if(strcmp1(&filab[2697],"ME  ")==0)
 259     for(l=0;l<6**nk;l++){emnt[l]=emn[l];};
 260   
 261   for(jj=0;jj<*mcs;jj++){
 262     is=cs[17*jj+4];
 263     for(i=1;i<is;i++){
 264     
 265       for(l=0;l<*nk;l++){inumt[l+i**nk]=inum[l];}
 266     
 267       if((strcmp1(&filab[0],"U ")==0)||
 268          ((strcmp1(&filab[87],"NT  ")==0)&&(*ithermal>=2))){
 269         for(l1=0;l1<*nk;l1++){
 270           if(inocs[l1]==jj){
 271             for(l2=0;l2<4;l2++){
 272               l=mt*l1+l2;
 273               vt[l+mt**nk*i]=v[l];
 274             }
 275           }
 276         }
 277       }
 278     
 279       if((strcmp1(&filab[87],"NT  ")==0)&&(*ithermal<2)){
 280         for(l=0;l<*nk;l++){
 281           if(inocs[l]==jj) t1t[l+*nk*i]=t1[l];
 282         }
 283       }
 284     
 285       if(strcmp1(&filab[174],"S   ")==0){
 286         for(l1=0;l1<*nk;l1++){
 287           if(inocs[l1]==jj){
 288             for(l2=0;l2<6;l2++){
 289               l=6*l1+l2;
 290               stnt[l+6**nk*i]=stn[l];
 291             }
 292           }
 293         }
 294       }
 295     
 296       if(strcmp1(&filab[261],"E   ")==0){
 297         for(l1=0;l1<*nk;l1++){
 298           if(inocs[l1]==jj){
 299             for(l2=0;l2<6;l2++){
 300               l=6*l1+l2;
 301               eent[l+6**nk*i]=een[l];
 302             }
 303           }
 304         }
 305       }
 306     
 307       if((strcmp1(&filab[348],"RF  ")==0)||(strcmp1(&filab[783],"RFL ")==0)){
 308         for(l1=0;l1<*nk;l1++){
 309           if(inocs[l1]==jj){
 310             for(l2=0;l2<4;l2++){
 311               l=mt*l1+l2;
 312               fnt[l+mt**nk*i]=fn[l];
 313             }
 314           }
 315         }
 316       }
 317     
 318       if(strcmp1(&filab[435],"PEEQ")==0){
 319         for(l=0;l<*nk;l++){
 320           if(inocs[l]==jj) epnt[l+*nk*i]=epn[l];
 321         }
 322       } 
 323     
 324       if(strcmp1(&filab[522],"ENER")==0){
 325         for(l=0;l<*nk;l++){
 326           if(inocs[l]==jj) enernt[l+*nk*i]=enern[l];
 327         }
 328       } 
 329     
 330       if(strcmp1(&filab[609],"SDV ")==0){
 331         for(l1=0;l1<*nk;l1++){
 332           if(inocs[l1]==jj){
 333             for(l2=0;l2<*nstate_;l2++){
 334               l=*nstate_*l1+l2;
 335               xstatent[l+*nstate_**nk*i]=xstaten[l];
 336             }
 337           } 
 338         }
 339       }
 340     
 341       if(strcmp1(&filab[696],"HFL ")==0){
 342         for(l1=0;l1<*nk;l1++){
 343           if(inocs[l1]==jj){
 344             for(l2=0;l2<3;l2++){
 345               l=3*l1+l2;
 346               qfnt[l+3**nk*i]=qfn[l];
 347             }
 348           }
 349         }
 350       }
 351     
 352       if(strcmp1(&filab[2697],"ME  ")==0){
 353         for(l1=0;l1<*nk;l1++){
 354           if(inocs[l1]==jj){
 355             for(l2=0;l2<6;l2++){
 356               l=6*l1+l2;
 357               emnt[l+6**nk*i]=emn[l];
 358             }
 359           }
 360         }
 361       }
 362     }
 363   }
 364   
 365   icntrl=-2;
 366   
 367   FORTRAN(rectcyl,(cot,vt,fnt,stnt,qfnt,eent,cs,&nkt,&icntrl,t,filab,
 368                    &imag,mi,emn));
 369   
 370   if(strcmp1(&filab[1044],"ZZS")==0){
 371       NNEW(neigh,ITG,40*net);
 372       NNEW(ipneigh,ITG,nkt);
 373   }
 374 
 375   frd(cot,&nkt,kont,ipkont,lakont,&net0,vt,stnt,inumt,nmethod,
 376             kode,filab,eent,t1t,fnt,time,epnt,ielmatt,matname,enernt,xstatent,
 377             nstate_,istep,iinc,ithermal,qfnt,&mode,noddiam,trab,inotrt,
 378             ntrans,orab,ielorien,norien,description,ipneigh,neigh,
 379             mi,stit,vr,vi,stnr,stni,vmax,stnmax,&ngraph,veold,ener,&net,
 380             cs,set,nset,istartset,iendset,ialset,eenmax,fnr,fni,emnt,
 381             thicke,jobnamec,output,qfx,cdn,mortar,cdnr,cdni,nmat);
 382 
 383   if(strcmp1(&filab[1044],"ZZS")==0){SFREE(ipneigh);SFREE(neigh);}
 384   
 385   if((strcmp1(&filab[0],"U ")==0)||
 386      ((strcmp1(&filab[87],"NT  ")==0)&&(*ithermal>=2))) SFREE(vt);
 387   if((strcmp1(&filab[87],"NT  ")==0)&&(*ithermal<2)) SFREE(t1t);
 388   if((strcmp1(&filab[174],"S   ")==0)||(strcmp1(&filab[1044],"ZZS ")==0)||
 389      (strcmp1(&filab[1044],"ERR ")==0)) 
 390      SFREE(stnt);
 391   if(strcmp1(&filab[261],"E   ")==0) SFREE(eent);
 392   if((strcmp1(&filab[348],"RF  ")==0)||(strcmp1(&filab[783],"RFL ")==0))
 393         SFREE(fnt);
 394   if(strcmp1(&filab[435],"PEEQ")==0) SFREE(epnt);
 395   if(strcmp1(&filab[522],"ENER")==0) SFREE(enernt);
 396   if(strcmp1(&filab[609],"SDV ")==0) SFREE(xstatent);
 397   if(strcmp1(&filab[696],"HFL ")==0) SFREE(qfnt);
 398   if((strcmp1(&filab[1044],"ZZS ")==0)||(strcmp1(&filab[1044],"ERR ")==0)||
 399      (strcmp1(&filab[2175],"CONT")==0)) SFREE(stit);
 400   if(strcmp1(&filab[2697],"ME  ")==0) SFREE(emnt);
 401 
 402   SFREE(kont);SFREE(ipkont);SFREE(lakont);SFREE(ielmatt);
 403   SFREE(inumt);SFREE(cot);if(*ntrans>0)SFREE(inotrt);
 404   SFREE(inocs);SFREE(ielcs);
 405   return;
 406 }
 407 

/* [<][>][^][v][top][bottom][index][help] */