~ubuntu-branches/ubuntu/precise/code-saturne/precise

« back to all changes in this revision

Viewing changes to salome/libmilieu/runmilieu.c

  • Committer: Package Import Robot
  • Author(s): Sylvestre Ledru
  • Date: 2011-11-01 17:43:32 UTC
  • mto: (6.1.7 sid)
  • mto: This revision was merged to the branch mainline in revision 11.
  • Revision ID: package-import@ubuntu.com-20111101174332-tl4vk45no0x3emc3
Tags: upstream-2.1.0
ImportĀ upstreamĀ versionĀ 2.1.0

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/*
 
2
  This file is part of Code_Saturne, a general-purpose CFD tool.
 
3
 
 
4
  Copyright (C) 1998-2011 EDF S.A.
 
5
 
 
6
  This program is free software; you can redistribute it and/or modify it under
 
7
  the terms of the GNU General Public License as published by the Free Software
 
8
  Foundation; either version 2 of the License, or (at your option) any later
 
9
  version.
 
10
 
 
11
  This program is distributed in the hope that it will be useful, but WITHOUT
 
12
  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
 
13
  FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
 
14
  details.
 
15
 
 
16
  You should have received a copy of the GNU General Public License along with
 
17
  this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
 
18
  Street, Fifth Floor, Boston, MA 02110-1301, USA.
 
19
*/
 
20
 
 
21
/*----------------------------------------------------------------------------*/
 
22
 
 
23
/*----------------------------------------------------------------------------
 
24
 * Standard C library headers
 
25
 *----------------------------------------------------------------------------*/
 
26
 
 
27
#include <stdio.h>
 
28
 
 
29
/*----------------------------------------------------------------------------
 
30
 * Local headers
 
31
 *----------------------------------------------------------------------------*/
 
32
 
 
33
#include "utilitaires.h"
 
34
#include "donnees.h"
 
35
#include "communication.h"
 
36
 
 
37
/*----------------------------------------------------------------------------
 
38
 *  Header for the current file
 
39
 *----------------------------------------------------------------------------*/
 
40
 
 
41
#include "runmilieu.h"
 
42
 
 
43
/*----------------------------------------------------------------------------*/
 
44
 
 
45
#ifdef __cplusplus
 
46
extern "C" {
 
47
#endif
 
48
 
 
49
/*============================================================================
 
50
 *  Global variables
 
51
 *============================================================================*/
 
52
 
 
53
int     nb_dyn = 0;
 
54
int     nb_for = 0;
 
55
int     ntcast = 0;
 
56
double  lref = 0.0;
 
57
 
 
58
double  *xast = NULL;
 
59
double  *xvast = NULL;
 
60
double  *xvasa = NULL;
 
61
double  *xastp = NULL;
 
62
 
 
63
double  *foras = NULL;
 
64
double  *foaas = NULL;
 
65
double  *fopas = NULL;
 
66
 
 
67
/*============================================================================
 
68
 * Public function definitions
 
69
 *============================================================================*/
 
70
 
 
71
/*----------------------------------------------------------------------------
 
72
 * "runmilieu" function
 
73
 *----------------------------------------------------------------------------*/
 
74
 
 
75
void
 
76
runmilieu(void *icompo)
 
77
{
 
78
  /* local variables */
 
79
 
 
80
  int i,j;
 
81
  int ierr = 0;
 
82
  int icv = 0;
 
83
  double c1, c2, c3;
 
84
  double alpha, beta;
 
85
  double dt_ast, dt_sat;
 
86
  double dtold = 0.;
 
87
  double dt = dtref;
 
88
 
 
89
  /* Input data for the coupled calculation are defined in SALOME by
 
90
     the study case's XML file */
 
91
 
 
92
  /* Initialize communication */
 
93
  if ((ierr = inicom(icompo)) >= 0) {
 
94
    printf(" Initializing communication\n");
 
95
  }
 
96
 
 
97
  /* Send parameters to codes */
 
98
  if ((ierr = send_param(icompo)) >= 0) {
 
99
    printf(" Send calculation parameters to codes\n");
 
100
  }
 
101
 
 
102
  /* Compute array sizes and initialize */
 
103
 
 
104
  /* Receive geometric data (nb_for and nb_dyn) */
 
105
  if (ierr >= 0) {
 
106
    ierr = recv_geom(icompo);
 
107
 
 
108
    printf("----------------------------------\n");
 
109
    printf(" Geometric parameters\n");
 
110
    printf("   number of coupled faces: %i\n", nb_for);
 
111
    printf("   number of coupled nodes: %i\n", nb_dyn);
 
112
    printf("   reference length (m): %4.2le\n", lref  );
 
113
    printf("----------------------------------\n");
 
114
 
 
115
    /* dynamics */
 
116
    alldyn();
 
117
 
 
118
    /* efforts  */
 
119
    allfor();
 
120
 
 
121
    /* Prediction coefficients */
 
122
    c1    = 0.;
 
123
    c2    = 0.;
 
124
    c3    = 0.;
 
125
    beta  = 0.;
 
126
    alpha = 0.;
 
127
  }
 
128
 
 
129
  /* Send geometric data to Code_Aster (remove in the future)*/
 
130
  if (ierr >= 0) {
 
131
    ierr = send_geom(icompo);
 
132
  }
 
133
 
 
134
  /* Initialize time step */
 
135
  dt = 0.;
 
136
  dt_ast = 0.;
 
137
  dt_sat = 0.;
 
138
 
 
139
  /* Initialize coupling iteration */
 
140
  ntcast = 0;
 
141
 
 
142
 /* Main loop */
 
143
  i = 1;
 
144
  while (ierr >= 0) {
 
145
    printf("\n");
 
146
    printf("----------------------------------------------------\n");
 
147
    printf("\n");
 
148
    printf("*********************************\n");
 
149
    printf("*         iteration %i          *\n", i);
 
150
    printf("*********************************\n");
 
151
    printf("\n");
 
152
    printf("\n");
 
153
 
 
154
    /* Info on time scheme */
 
155
    if (nbssit <= 1) {
 
156
      printf("Explicit time-stepping scheme\n");
 
157
    }
 
158
    else {
 
159
      printf("Implicit time-stepping scheme\n");
 
160
      printf("  number of sub-iterations: %i\n", nbssit);
 
161
    }
 
162
 
 
163
    /* Manage time steps */
 
164
 
 
165
    /* Receive time steps from Code_Aster and Code_Saturne */
 
166
    ierr = recv_pdt(icompo,&(dt_ast), &(dt_sat), i);
 
167
 
 
168
    printf("----------------------------------\n");
 
169
    printf("reference time step:     %4.21e \n", dtref );
 
170
    printf("Code_Saturne time step:  %4.2le \n", dt_sat);
 
171
    printf("Code_Aster time step:    %4.2le \n", dt_ast);
 
172
 
 
173
    /* choose the smallest time step: dt = dt_ast; */
 
174
    dt = dtref;
 
175
    if (dt > dt_ast) {
 
176
      dt = dt_ast;
 
177
    }
 
178
    if (dt > dt_sat) {
 
179
      dt = dt_sat;
 
180
    }
 
181
 
 
182
    /* Send the selected time step */
 
183
    if (ierr >= 0) ierr = send_pdt(icompo, dt, i);
 
184
 
 
185
    printf("selected time step:      %4.2le \n", dt);
 
186
    printf("----------------------------------\n");
 
187
    printf("\n\n");
 
188
 
 
189
    icv = 0;
 
190
 
 
191
    j = 1;
 
192
    while (ierr >= 0) {
 
193
 
 
194
      printf("*********************************\n");
 
195
      printf("*     sub - iteration %i        *\n", j);
 
196
      printf("*********************************\n");
 
197
      printf("\n\n");
 
198
 
 
199
      /* increment coupling iteration */
 
200
      printf("midde\n");
 
201
      ntcast = ntcast + 1;
 
202
      printf("ntcast = %i\n", ntcast);
 
203
 
 
204
      /* printf("***************************************\n"); */
 
205
      /* printf("*        predict displacements        *\n"); */
 
206
      /* printf("***************************************\n"); */
 
207
 
 
208
      /* Predict displacements */
 
209
 
 
210
      c1 = 0.;
 
211
      c2 = 0.;
 
212
      c3 = 0.;
 
213
 
 
214
      /* seperate prediction for explicit/implicit cases */
 
215
      if (j == 1) {
 
216
        alpha = 0.5;
 
217
        beta  = 0.;
 
218
        c1    = 1.;
 
219
        c2    = (alpha + beta) * dt ;
 
220
        c3    = -beta * dtold ;
 
221
        pred(xastp, xast, xvast, xvasa, c1, c2, c3, nb_dyn);
 
222
      }
 
223
 
 
224
      if (j > 1) {
 
225
        alpha = 0.5;
 
226
        c1    = alpha;
 
227
        c2    = 1. - alpha ;
 
228
        c3    = 0.;
 
229
        pred(xastp, xast, xastp, xast, c1, c2, c3, nb_dyn);
 
230
      }
 
231
 
 
232
      printf("--------------------------------------------\n");
 
233
      printf("Displacement prediction coefficients\n");
 
234
      printf(" C1: %4.2le\n", c1);
 
235
      printf(" C2: %4.2le\n", c2);
 
236
      printf(" C3: %4.2le\n", c3);
 
237
      printf("--------------------------------------------\n");
 
238
      printf("\n\n");
 
239
 
 
240
      /* send predicted displacements */
 
241
      if (ierr >= 0) ierr = send_dyn(icompo);
 
242
 
 
243
      /* explicit case: no need for a convergence test */
 
244
 
 
245
      /* implicit case: needs a convergence test */
 
246
 
 
247
      /* printf("***************************************\n"); */
 
248
      /* printf("*  end of displacements prediction    *\n"); */
 
249
      /* printf("***************************************\n"); */
 
250
 
 
251
      /* printf("*********************************\n"); */
 
252
      /* printf("*       forces prediction       *\n"); */
 
253
      /* printf("*********************************\n"); */
 
254
 
 
255
      /* Receive forces */
 
256
      ierr = recv_for(icompo);
 
257
 
 
258
      /* No difference between explicit and implicit cases for forces */
 
259
      alpha = 2.0;
 
260
      c1    = alpha;
 
261
      c2    = 1-alpha;
 
262
      c3    = 0.;
 
263
      pred(fopas, foras, foaas, foaas, c1, c2, c3, nb_for);
 
264
      printf("--------------------------------------\n");
 
265
      printf("Forces prediction coefficients\n");
 
266
      printf(" C1: %4.2le\n",c1);
 
267
      printf(" C2: %4.2le\n",c2);
 
268
      printf(" C3: %4.2le\n",c3);
 
269
      printf("--------------------------------------\n");
 
270
      printf("\n\n");
 
271
 
 
272
      /* send des forces */
 
273
      if (ierr >= 0) ierr = send_for(icompo);
 
274
 
 
275
      /* printf("*********************************\n"); */
 
276
      /* printf("*   end of forces prediction    *\n"); */
 
277
      /* printf("*********************************\n"); */
 
278
 
 
279
      printf("\n");
 
280
 
 
281
      /* explicit case: no need fo a convergence test */
 
282
      if (nbssit <= 1) {
 
283
        /* handle convergence even when no test is done */
 
284
        icv =1;
 
285
        if (ierr >= 0) ierr = send_icv1(icompo,icv);
 
286
        if (ierr >= 0) ierr = recv_icv(icompo,&(icv));
 
287
        icv = 1;
 
288
        if (ierr >= 0) ierr = send_icv2(icompo,icv);
 
289
 
 
290
        /* receive displacements effectively calculated by Code_Aster */
 
291
        if (ierr >= 0) ierr = recv_dyn(icompo);
 
292
 
 
293
        /* record previous values */
 
294
        val_ant();
 
295
 
 
296
        break;
 
297
      }
 
298
 
 
299
      /* implicit case: requires a convergence test */
 
300
      else {
 
301
        /* compute icv */
 
302
        if (ierr >= 0) ierr = conv(&(icv));
 
303
        if (ierr >= 0) ierr = send_icv1(icompo,icv);
 
304
        if (ierr >= 0) ierr = recv_icv(icompo,&(icv));
 
305
        if (ierr >= 0) ierr = send_icv2(icompo,icv);
 
306
 
 
307
        if((j>=nbssit) || (icv == 1)) {
 
308
          /* receive displacements effectivemey computed by Code_Aster */
 
309
          /* Receive displacements */
 
310
          if (ierr >= 0) ierr = recv_dyn(icompo);
 
311
 
 
312
          /* then send to Code_Saturne ? the question remains open... */
 
313
          /* if necessary, function to send these displs. should be
 
314
             created in middle and matching receive in Code_Saturne */
 
315
          /* if (ierr >= 0) ierr = send2_dyn(); */
 
316
 
 
317
          /* receive displacements effectiveley calculated by Code_Aster */
 
318
          if (ierr >= 0) ierr = recv_dyn(icompo);
 
319
          break;
 
320
        }
 
321
        else {
 
322
          j = j+1;
 
323
        }
 
324
      }
 
325
    } /* end of sub-iterations loop */
 
326
 
 
327
    /* iterations test */
 
328
    if (i >= nbpdtm) {
 
329
      ierr = -1;
 
330
    }
 
331
    /* end of iterations test */
 
332
 
 
333
    i = i+1;
 
334
 
 
335
    /* save time step */
 
336
    dtold = dt;
 
337
 
 
338
  } /* en of iterations loop */
 
339
 
 
340
  ierr = calfin(icompo);
 
341
}
 
342
 
 
343
/*----------------------------------------------------------------------------*/
 
344
 
 
345
#ifdef __cplusplus
 
346
}
 
347
#endif