00001 #include "cpnpsolvert.h"
00002 #include "globalt.h"
00003 #include "globmatt.h"
00004 #include "transprint.h"
00005 #include <string.h>
00006
00007
00008
00009
00010
00011
00012
00013 void solve_nonstationary_growing_problem ()
00014 {
00015 long i,j,n,nce,ncd,lcid;
00016 double s,dt,prev_time,end_time,alpha,*d,*p,*lhs,*tdlhs,*rhs;
00017
00018 n=Ndoft;
00019 lcid=0;
00020
00021 lhs = Lsrst->give_lhs (0);
00022 tdlhs = Lsrst->give_tdlhs (0);
00023 rhs = Lsrst->give_rhs (0);
00024
00025 d = new double [n];
00026 p = new double [n];
00027
00028
00029 nullv (lhs,n);
00030 nullv (tdlhs,n);
00031 nullv (d,n);
00032 nullv (p,n);
00033
00034 alpha=Tp->alpha;
00035
00036
00037 Tp->time=Tp->timecont.starttime ();
00038 prev_time=Tp->time;
00039
00040 end_time = Tp->timecont.endtime ();
00041
00042 dt = Tp->timecont.initialtimeincr ();
00043
00044
00045
00046
00047 i=0;
00048
00049
00050 approximation ();
00051 Tm->initmaterialmodels();
00052 print_initt(-1, "wt");
00053 print_stept(0,i,Tp->time,NULL);
00054
00055
00056 Tt->lhs_save (lhs,Lsrst->lhsi,tdlhs);
00057
00058 do{
00059
00060 if (Mesprt != 0) fprintf (stdout,"\n\n increment number %ld, time %f, time step = %f\n",i,Tp->time,dt);
00061
00062
00063
00064
00065 nce = Gtt->search_newelem (Tp->time,prev_time);
00066
00067
00068 Gtt->update_elem (Tp->time);
00069
00070
00071
00072
00073 ncd = Gtt->search_newdofs (Tp->time,prev_time);
00074
00075
00076 Gtt->update_dofs (Tp->time);
00077
00078
00079 Ndoft = Gtt->codenum_generation (Outt);
00080 n=Ndoft;
00081
00082 Tm->updateipval ();
00083 i++;
00084
00085 if (nce!=0 || ncd!=0){
00086
00087 Tt->lhs_restore (lhs,Lsrst->lhsi,tdlhs);
00088 }
00089
00090
00091 if (Kmat != NULL){
00092 delete Kmat;
00093 Kmat=NULL;
00094 }
00095
00096 if (Cmat != NULL){
00097 delete Cmat;
00098 Cmat=NULL;
00099 }
00100
00101
00102
00103
00104 if (nce>0){
00105 Tt->initial_nodval ();
00106 }
00107
00108
00109 capacity_matrix (0);
00110
00111
00112 conductivity_matrix (0);
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123 for (j=0;j<n;j++){
00124 p[j]=lhs[j] + (1.0-alpha)*dt*tdlhs[j];
00125 }
00126
00127
00128 Kmat->gmxv (p,d);
00129
00130 Kmat->scalgm (dt*alpha);
00131
00132 Kmat->addgm (1.0,*Cmat);
00133
00134
00135
00136 trfel_right_hand_side (0,rhs,n);
00137
00138 for (j=0;j<n;j++){
00139 rhs[j] = rhs[j] - d[j];
00140 d[j]=tdlhs[j];
00141 }
00142
00143
00144 Tp->ssle->solve_system (Gtt,Kmat,tdlhs,rhs,Outt);
00145
00146 for (j=0;j<n;j++){
00147 s=(1.0-alpha)*d[j]+alpha*tdlhs[j];
00148 lhs[j]+=dt*s;
00149 }
00150
00151
00152
00153 solution_correction ();
00154
00155 Tt->lhs_save (lhs,Lsrst->lhsi,tdlhs);
00156
00157 approximation ();
00158
00159
00160 prev_time=Tp->time;
00161 Tp->time=Tp->timecont.newtime ();
00162 dt = Tp->timecont.actualbacktimeincr ();
00163
00164 print_stept(0,i,Tp->time,NULL);
00165 print_flusht();
00166
00167 }while(Tp->time<end_time);
00168
00169 delete [] p; delete [] d;
00170 print_closet();
00171 }
00172
00173
00174
00175
00176
00177
00178
00179
00180 void solve_nonstationary_growing_vform ()
00181 {
00182 long i,j,n,nce,ncd,lcid;
00183 double dt,prev_time,end_time,alpha,*d,*p,*lhs,*tdlhs,*rhs;
00184
00185 n=Ndoft;
00186 lcid=0;
00187
00188 lhs = Lsrst->give_lhs (0);
00189 tdlhs = Lsrst->give_tdlhs (0);
00190 rhs = Lsrst->give_rhs (0);
00191
00192 d = new double [n];
00193 p = new double [n];
00194
00195
00196 nullv (lhs,n);
00197 nullv (tdlhs,n);
00198 nullv (d,n);
00199 nullv (p,n);
00200
00201 alpha=Tp->alpha;
00202
00203
00204 Tp->time=Tp->timecont.starttime ();
00205 prev_time=Tp->time;
00206
00207 end_time = Tp->timecont.endtime ();
00208
00209 dt = Tp->timecont.initialtimeincr ();
00210
00211
00212
00213
00214 i=0;
00215
00216
00217 approximation ();
00218 Tm->initmaterialmodels();
00219 print_initt(-1, "wt");
00220 print_stept(0,i,Tp->time,NULL);
00221
00222
00223 Tt->lhs_save (lhs,Lsrst->lhsi,tdlhs);
00224
00225 do{
00226
00227 if (Mesprt != 0) fprintf (stdout,"\n\n increment number %ld, time %f, time step = %f\n",i,Tp->time,dt);
00228
00229
00230
00231
00232 nce = Gtt->search_newelem (Tp->time,prev_time);
00233
00234
00235 Gtt->update_elem (Tp->time);
00236
00237
00238
00239
00240 ncd = Gtt->search_newdofs (Tp->time,prev_time);
00241
00242
00243 Gtt->update_dofs (Tp->time);
00244
00245
00246 Ndoft = Gtt->codenum_generation (Outt);
00247 n=Ndoft;
00248
00249 Tm->updateipval ();
00250 i++;
00251
00252 if (nce!=0 || ncd!=0){
00253
00254 Tt->lhs_restore (lhs,Lsrst->lhsi,tdlhs);
00255 }
00256
00257
00258 if (Kmat != NULL){
00259 delete Kmat;
00260 Kmat=NULL;
00261 }
00262
00263 if (Cmat != NULL){
00264 delete Cmat;
00265 Cmat=NULL;
00266 }
00267
00268
00269
00270
00271 if (nce>0){
00272 Tt->initial_nodval ();
00273 }
00274
00275
00276 capacity_matrix (0);
00277
00278
00279 conductivity_matrix (0);
00280
00281
00282
00283
00284
00285
00286
00287
00288
00289
00290 for (j=0;j<n;j++){
00291 d[j]=lhs[j] + (1.0-alpha)*dt*tdlhs[j];
00292 }
00293
00294
00295 Kmat->gmxv (d,p);
00296
00297 Kmat->scalgm (dt*alpha);
00298
00299 Kmat->addgm (1.0,*Cmat);
00300
00301
00302
00303 trfel_right_hand_side (0,rhs,n);
00304
00305 for (j=0;j<n;j++){
00306 rhs[j] = rhs[j] - p[j];
00307 }
00308
00309
00310 Tp->ssle->solve_system (Gtt,Kmat,tdlhs,rhs,Outt);
00311
00312 for (j=0;j<n;j++){
00313 lhs[j] = d[j] + alpha*dt*tdlhs[j];
00314 }
00315
00316
00317 solution_correction ();
00318
00319 Tt->lhs_save (lhs,Lsrst->lhsi,tdlhs);
00320
00321 approximation ();
00322
00323
00324 prev_time=Tp->time;
00325 Tp->time=Tp->timecont.newtime ();
00326 dt = Tp->timecont.actualbacktimeincr ();
00327
00328 print_stept(0,i,Tp->time,NULL);
00329 print_flusht();
00330
00331 }while(Tp->time<end_time);
00332
00333 delete [] p; delete [] d;
00334 print_closet();
00335 }
00336
00337