APBS  3.0.0
mgfasd.c
1 
55 #include "mgfasd.h"
56 
57 VPUBLIC void Vfmvfas(int *nx, int *ny, int *nz,
58  double *x,
59  int *iz,
60  double *w0, double *w1, double *w2, double *w3, double *w4,
61  int *istop, int *itmax, int *iters, int *ierror,
62  int *nlev, int *ilev, int *nlev_real,
63  int *mgsolv, int *iok, int *iinfo,
64  double *epsiln, double *errtol, double *omega,
65  int *nu1, int *nu2, int *mgsmoo,
66  int *ipc, double *rpc,
67  double *pc, double *ac, double *cc, double *fc, double *tru) {
68 
69  // Other Declarations
70  int level, itmxd, nlevd, iterd, iokd;
71  int nxf, nyf, nzf;
72  int nxc, nyc, nzc;
73  int istpd, iinfod;
74  double errd;
75 
76  // Utility variables
77  int numlev;
78 
79  MAT2(iz, 50, *nlev);
80 
81  // Recover gridsizes
82  nxf = *nx;
83  nyf = *ny;
84  nzf = *nz;
85 
86  numlev = *nlev - 1;
87  Vmkcors(&numlev, &nxf, &nyf, &nzf, &nxc, &nyc, &nzc);
88 
89  // Move up grids: interpolate solution to finer, do v cycle
90  if (*iinfo != 0) {
91 
92  Vnm_print(2, "Vfmvfas: starting: (%d,%d,%d) (%d,%d,%d)\n",
93  nxf, nyf, nzf, nxc, nyc, nzc);
94  }
95 
96  for (level = *nlev_real; level >= *ilev + 1; level--) {
97 
98  // Call mv cycle
99  errd = 1.0e-5;
100  itmxd = 1;
101  nlevd = *nlev_real - level + 1;
102  iterd = 0;
103  iokd = 2;
104  iinfod = *iinfo;
105  istpd = *istop;
106  if (*iinfo >= 2)
107  iokd = 2;
108 
109  Vmvfas(&nxc, &nyc, &nzc,
110  x,
111  iz,
112  w0, w1, w2, w3, w4,
113  &istpd, &itmxd, &iterd, ierror,
114  &nlevd, &level, nlev_real,
115  mgsolv, &iokd, &iinfod,
116  epsiln, errtol, omega,
117  nu1, nu2, mgsmoo,
118  ipc, rpc,
119  pc, ac, cc, fc, tru);
120 
121  // Find new grid size
122  numlev = 1;
123  Vmkfine(&numlev, &nxc, &nyc, &nzc, &nxf, &nyf, &nzf);
124 
125  // Interpolate to next finer grid
126  VinterpPMG(&nxc, &nyc, &nzc,
127  &nxf, &nyf, &nzf,
128  RAT( x, VAT2(iz, 1, level)),
129  RAT( x, VAT2(iz, 1, level-1)),
130  RAT(pc, VAT2(iz, 11, level-1)));
131 
132  // New grid size
133  nxc = nxf;
134  nyc = nyf;
135  nzc = nzf;
136  }
137 
138 
139 
140  // Call mv cycle
141  level = *ilev;
142 
143  Vmvfas(&nxf, &nyf, &nzf,
144  x, iz,
145  w0, w1, w2, w3, w4,
146  istop, itmax, iters,
147  ierror, nlev, &level, nlev_real,
148  mgsolv, iok, iinfo,
149  epsiln, errtol, omega,
150  nu1, nu2, mgsmoo,
151  ipc, rpc,
152  pc, ac, cc, fc, tru);
153 }
154 
155 
156 
157 VPUBLIC void Vmvfas(int *nx, int *ny, int *nz,
158  double *x,
159  int *iz,
160  double *w0, double *w1, double *w2, double *w3, double *w4,
161  int *istop, int *itmax, int *iters, int *ierror,
162  int *nlev, int *ilev, int *nlev_real,
163  int *mgsolv, int *iok, int *iinfo,
164  double *epsiln, double *errtol, double *omega,
165  int *nu1, int *nu2, int *mgsmoo,
166  int *ipc, double *rpc,
167  double *pc, double *ac, double *cc, double *fc, double *tru) {
168 
169  // Other declarations
170  int level, lev;
171  int itmax_s, iters_s, nuuu, ivariv, mgsmoo_s, iresid;
172  int nxf, nyf, nzf;
173  int nxc, nyc, nzc;
174  int iadjoint;
175  double errtol_s;
176  double rsden, rsnrm, orsnrm;
177  double xdamp;
178 
179  int numlev;
180  double alpha;
181 
182  MAT2(iz, 50, *nlev);
183 
184  WARN_UNTESTED;
185 
186  // Recover level information
187  level = 1;
188  lev = (*ilev - 1) + level;
189 
190  // Recover gridsizes
191  nxf = *nx;
192  nyf = *ny;
193  nzf = *nz;
194 
195  numlev = *nlev - 1;
196  Vmkcors(&numlev, &nxf, &nyf, &nzf, &nxc, &nyc, &nzc);
197 
198  // Do some i/o if requested
199  if (*iinfo != 0) {
200  Vnm_print(2, "Vmvfas: starting: (%d, %d, %d) (%d, %d, %d)\n",
201  nxf, nyf, nzf, nxc, nyc, nzc);
202  }
203 
204  // Initial wall clock
205  if (*iok != 0) {
206  Vprtstp(*iok, -1, 0.0, 0.0, 0.0);
207  }
208 
209  /**************************************************************
210  *** note: if (iok != 0) then: use a stopping test. ***
211  *** else: use just the itmax to stop iteration. ***
212  **************************************************************
213  *** istop=0 most efficient (whatever it is) ***
214  *** istop=1 relative residual ***
215  *** istop=2 rms difference of successive iterates ***
216  *** istop=3 relative true error (provided for testing) ***
217  **************************************************************/
218 
219  // Compute denominator for stopping criterion
220  if (*iok != 0) {
221  if (*istop == 0) {
222  rsden = 1.0;
223  } else if (*istop == 1) {
224 
225  // Compute initial residual with zero initial guess
226  // this is analogous to the linear case where one can
227  // simply take norm of rhs for a zero initial guess
228  Vazeros(&nxf, &nyf, &nzf, w1);
229 
230  Vnmresid(&nxf, &nyf, &nzf,
231  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
232  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
233  RAT( fc, VAT2(iz, 1, lev)),
234  w1, w2, w3);
235 
236  rsden = Vxnrm1(&nxf, &nyf, &nzf, w2);
237 
238  } else if (*istop == 2) {
239  rsden = VSQRT(nxf * nyf * nzf);
240  } else if (*istop == 3) {
241  rsden = Vxnrm2(&nxf, &nyf, &nzf, RAT(tru, VAT2(iz, 1, lev)));
242  } else if (*istop == 4) {
243  rsden = Vxnrm2(&nxf, &nyf, &nzf, RAT(tru, VAT2(iz, 1, lev)));
244  } else if (*istop == 5) {
245  Vnmatvec(&nxf, &nyf, &nzf,
246  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
247  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
248  RAT(tru, VAT2(iz, 1, lev)),
249  w1, w2);
250  rsden = VSQRT(Vxdot(&nxf, &nyf, &nzf,RAT(tru, VAT2(iz, 1, lev)), w1));
251  } else {
252  Vnm_print(2, "Vmvfas: bad istop value: %d\n", *istop);
253  }
254  if (rsden == 0.0) {
255  rsden = 1.0;
256  Vnm_print(2, "Vmfas: rhs is zero on finest level\n");
257  }
258  rsnrm = rsden;
259  orsnrm = rsnrm;
260 
261  Vprtstp(*iok, 0, rsnrm, rsden, orsnrm);
262  }
263 
264 
265 
266  /********************************************************************
267  *** solve directly if nlev = 1
268  ********************************************************************/
269 
270  // Solve directly if on the coarse grid
271  if (*nlev == 1) {
272 
273  //solve with ncghs, mgsmoo_s=4 (no residual)
274  iresid = 0;
275  iadjoint = 0;
276  itmax_s = 100;
277  iters_s = 0;
278  errtol_s = *epsiln;
279  mgsmoo_s = *mgsmoo;
280  Vazeros(&nxf, &nyf, &nzf,RAT(x, VAT2(iz, 1, lev)));
281  Vnsmooth(&nxf, &nyf, &nzf,
282  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
283  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
284  RAT( fc, VAT2(iz, 1, lev)),
285  RAT( x, VAT2(iz, 1, lev)),
286  w1, w2, w3,
287  &itmax_s, &iters_s, &errtol_s, omega,
288  &iresid, &iadjoint, &mgsmoo_s);
289 
290  // Compute the stopping test
291  *iters = 1;
292  if (*iok != 0) {
293  orsnrm = rsnrm;
294  if (*istop == 0) {
295  Vnmresid(&nxf, &nyf, &nzf,
296  RAT(ipc, VAT2(iz, 5, lev)),
297  RAT(rpc, VAT2(iz, 6, lev)),
298  RAT( ac, VAT2(iz, 7, lev)),
299  RAT( cc, VAT2(iz, 1, lev)),
300  RAT( fc, VAT2(iz, 1, lev)),
301  RAT( x, VAT2(iz, 1, lev)),
302  w1, w2);
303  rsnrm = Vxnrm1(&nxf, &nyf, &nzf, w1);
304  } else if (*istop == 1) {
305  Vnmresid(&nxf, &nyf, &nzf,
306  RAT(ipc, VAT2(iz, 5, lev)),
307  RAT(rpc, VAT2(iz, 6, lev)),
308  RAT( ac, VAT2(iz, 7, lev)),
309  RAT( cc, VAT2(iz, 1, lev)),
310  RAT( fc, VAT2(iz, 1, lev)),
311  RAT( x, VAT2(iz, 1, lev)),
312  w1, w2);
313  rsnrm = Vxnrm1(&nxf, &nyf, &nzf, w1);
314  } else if (*istop == 2) {
315  Vxcopy(&nxf, &nyf, &nzf,
316  RAT(tru, VAT2(iz, 1, lev)),
317  w1);
318  alpha = -1.0;
319  Vxaxpy(&nxf, &nyf, &nzf,
320  &alpha, RAT(x, VAT2(iz, 1, lev)), w1);
321  rsnrm = Vxnrm1(&nxf, &nyf, &nzf, w1);
322  Vxcopy(&nxf, &nyf, &nzf,
323  RAT(x, VAT2(iz, 1, lev)),
324  RAT(tru, VAT2(iz, 1, lev)));
325  } else if (*istop == 3) {
326  Vxcopy(&nxf, &nyf, &nzf,
327  RAT(tru, VAT2(iz, 1, lev)),
328  w1);
329  alpha = -1.0;
330  Vxaxpy(&nxf, &nyf, &nzf,
331  &alpha, RAT(x, VAT2(iz, 1, lev)), w1);
332  rsnrm = Vxnrm2(&nxf, &nyf, &nzf, w1);
333  } else if (*istop == 4) {
334  Vxcopy(&nxf, &nyf, &nzf,
335  RAT(tru, VAT2(iz, 1, lev)),
336  w1);
337  alpha = -1.0;
338  Vxaxpy(&nxf, &nyf, &nzf,
339  &alpha, RAT(x, VAT2(iz, 1, lev)), w1);
340  rsnrm = Vxnrm2(&nxf, &nyf, &nzf, w1);
341  } else if (*istop == 5) {
342  Vxcopy(&nxf, &nyf, &nzf,
343  RAT(tru, VAT2(iz, 1, lev)),
344  w1);
345  alpha = -1.0;
346  Vxaxpy(&nxf, &nyf, &nzf,
347  &alpha, RAT(x, VAT2(iz, 1, lev)), w1);
348  Vnmatvec(&nxf, &nyf, &nzf,
349  RAT(ipc, VAT2(iz, 5, lev)),
350  RAT(rpc, VAT2(iz, 6, lev)),
351  RAT( ac, VAT2(iz, 7, lev)),
352  RAT( cc, VAT2(iz, 1, lev)),
353  w1, w2, w3);
354  rsnrm = VSQRT(Vxdot(&nxf, &nyf, &nzf, w1, w2));
355  } else {
356  Vnm_print(2, "Vmvcs: bad istop value: %d\n", *istop);
357  }
358 
359  Vprtstp(*iok, *iters, rsnrm, rsden, orsnrm);
360  }
361  return;
362  }
363 
364  /*********************************************************************
365  *** begin mg iteration (note nxf,nyf,nzf changes during loop)
366  *********************************************************************/
367 
368  // setup for the v-cycle looping ***
369  *iters = 0;
370  while(1) {
371 
372  // Finest level initialization
373  level = 1;
374  lev = (*ilev - 1) + level;
375 
376  // Nu1 pre-smoothings on fine grid (with residual)
377  iresid = 1;
378  iadjoint = 0;
379  iters_s = 0;
380  errtol_s = 0.0;
381  nuuu = Vivariv(nu1, &lev);
382  Vnsmooth(&nxf, &nyf, &nzf,
383  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
384  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
385  RAT( fc, VAT2(iz, 1, lev)), RAT( x, VAT2(iz, 1, lev)),
386  w2, w3, w1,
387  &nuuu, &iters_s, &errtol_s, omega,
388  &iresid, &iadjoint, mgsmoo);
389  Vxcopy(&nxf, &nyf, &nzf, w1, RAT(w0, VAT2(iz, 1, lev)));
390 
391  /* *********************************************************************
392  * begin cycling down to coarse grid
393  * ********************************************************************/
394 
395  // Go down grids: restrict resid to coarser and smooth
396  for (level = 2; level <= *nlev; level++) {
397  lev = (*ilev - 1) + level;
398 
399  // Find new grid size
400  numlev = 1;
401  Vmkcors(&numlev, &nxf, &nyf, &nzf, &nxc, &nyc, &nzc);
402 
403  // Restrict residual to coarser grid
404  Vrestrc(&nxf, &nyf, &nzf, &nxc, &nyc, &nzc,
405  w1, RAT(w0, VAT2(iz, 1, lev)),
406  RAT(pc, VAT2(iz, 11, lev-1)));
407 
408  // Restrict (extract) solution to coarser grid
409  Vextrac(&nxf, &nyf, &nzf, &nxc, &nyc, &nzc,
410  RAT(x, VAT2(iz, 1, lev-1)), RAT(w4, VAT2(iz, 1, lev)));
411 
412  // New grid size
413  nxf = nxc;
414  nyf = nyc;
415  nzf = nzc;
416 
417  // Apply coarse grid operator to coarse grid soln
418  Vnmatvec(&nxf, &nyf, &nzf,
419  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
420  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
421  RAT( w4, VAT2(iz, 1, lev)), RAT( fc, VAT2(iz, 1, lev)),
422  w3);
423 
424  // Build coarse grid right hand side
425  alpha = 1.0;
426  Vxaxpy(&nxf, &nyf, &nzf, &alpha,
427  RAT(w0, VAT2(iz, 1, lev)), RAT(fc, VAT2(iz, 1, lev)));
428 
429  // If not on coarsest level yet...
430  if (level != *nlev) {
431 
432  // nu1 pre-smoothings on this level (with residual)
433  Vxcopy(&nxf, &nyf, &nzf,
434  RAT(w4, VAT2(iz, 1, lev)), RAT(x, VAT2(iz, 1, lev)));
435  iresid = 1;
436  iadjoint = 0;
437  iters_s = 0;
438  errtol_s = 0.0;
439  nuuu = Vivariv (nu1,&lev);
440  Vnsmooth(&nxf, &nyf, &nzf,
441  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
442  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
443  RAT( fc, VAT2(iz, 1, lev)),
444  RAT( x, VAT2(iz, 1, lev)),
445  w2, w3, w1,
446  &nuuu, &iters_s, &errtol_s, omega,
447  &iresid, &iadjoint, mgsmoo);
448  }
449 
450  // End of cycling down to coarse grid loop
451  }
452 
453 
454 
455  /* *********************************************************************
456  * begin coarse grid
457  * *********************************************************************/
458 
459  // Coarsest level
460  level = *nlev;
461  lev = (*ilev - 1) + level;
462 
463  // Solve on coarsest grid with ncghs, mgsmoo_s=4 (no residual)
464  iresid = 0;
465  iadjoint = 0;
466  itmax_s = 100;
467  iters_s = 0;
468  errtol_s = *epsiln;
469  mgsmoo_s = *mgsmoo;
470  Vazeros(&nxf, &nyf, &nzf, RAT(x, VAT2(iz, 1, lev)));
471  Vnsmooth(&nxf, &nyf, &nzf,
472  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
473  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
474  RAT( fc, VAT2(iz, 1, lev)),
475  RAT( x, VAT2(iz, 1, lev)),
476  w1, w2, w3,
477  &itmax_s, &iters_s, &errtol_s, omega,
478  &iresid, &iadjoint, &mgsmoo_s);
479 
480  /* *********************************************************************
481  * begin cycling back to fine grid
482  * ********************************************************************/
483 
484  // Move up grids: interpolate resid to finer and smooth
485  for (level = *nlev - 1; level >= 1; level--) {
486  lev = (*ilev - 1) + level;
487 
488  // Find new grid size
489  numlev = 1;
490  Vmkfine(&numlev, &nxf, &nyf, &nzf, &nxc, &nyc, &nzc);
491 
492  // Form difference of new approx at the coarse grid
493  alpha = -1.0;
494  Vxaxpy(&nxf, &nyf, &nzf, &alpha,
495  RAT(w4, VAT2(iz, 1, lev + 1)), RAT(x, VAT2(iz, 1, lev + 1)));
496 
497  // Call the line search (on the coarser level)
498  Vlinesearch(&nxf, &nyf, &nzf, &xdamp,
499  RAT(ipc, VAT2(iz, 5, lev + 1)),
500  RAT(rpc, VAT2(iz, 6, lev + 1)),
501  RAT( ac, VAT2(iz, 7, lev + 1)),
502  RAT( cc, VAT2(iz, 1, lev + 1)),
503  RAT( fc, VAT2(iz, 1, lev + 1)),
504  RAT( x, VAT2(iz, 1, lev + 1)),
505  RAT( w4, VAT2(iz, 1, lev + 1)),
506  RAT( w0, VAT2(iz, 1, lev + 1)),
507  w1, w2, w3);
508 
509  // Interpolate to next finer grid
510  VinterpPMG(&nxf, &nyf, &nzf, &nxc, &nyc, &nzc,
511  RAT(x, VAT2(iz, 1, lev + 1)),
512  w1,
513  RAT(pc, VAT2(iz, 11, lev)));
514 
515  // New grid size
516  nxf = nxc;
517  nyf = nyc;
518  nzf = nzc;
519 
520  // Perform the coarse grid correction
521  Vxaxpy(&nxf, &nyf, &nzf, &xdamp, w1, RAT(x, VAT2(iz, 1, lev)));
522 
523  // nu2 post-smoothings for correction (no residual) ***
524  iresid = 0;
525  iadjoint = 1;
526  iters_s = 0;
527  errtol_s = 0.0;
528  nuuu = Vivariv(nu2, &lev);
529  Vnsmooth(&nxf, &nyf, &nzf,
530  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
531  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
532  RAT( fc, VAT2(iz, 1, lev)),
533  RAT( x, VAT2(iz, 1, lev)),
534  w1, w2, w3,
535  &nuuu, &iters_s, &errtol_s, omega,
536  &iresid, &iadjoint, mgsmoo);
537  }
538 
539 
540  /* *********************************************************************
541  * iteration complete: do some i/o
542  * ********************************************************************/
543 
544  // Increment the iteration counter
545  iters = iters + 1;
546 
547  // Compute/check the current stopping test
548  if (*iok != 0) {
549  orsnrm = rsnrm;
550  if (*istop == 0) {
551  Vnmresid(&nxf, &nyf, &nzf,
552  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
553  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
554  RAT( fc, VAT2(iz, 1, lev)), RAT( x, VAT2(iz, 1, lev)),
555  w1, w2);
556  rsnrm = Vxnrm1(&nxf, &nyf, &nzf, w1);
557  } else if (*istop == 1) {
558  Vnmresid(&nxf, &nyf, &nzf,
559  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
560  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
561  RAT( fc, VAT2(iz, 1, lev)), RAT( x, VAT2(iz, 1, lev)),
562  w1, w2);
563  rsnrm = Vxnrm1(&nxf, &nyf, &nzf,w1);
564  } else if (*istop == 2) {
565  Vxcopy(&nxf, &nyf, &nzf, RAT(tru, VAT2(iz, 1, lev)), w1);
566  alpha = -1.0;
567  Vxaxpy(&nxf, &nyf, &nzf, &alpha, RAT(x, VAT2(iz, 1, lev)), w1);
568  rsnrm = Vxnrm1(&nxf, &nyf, &nzf, w1);
569  Vxcopy(&nxf, &nyf, &nzf,
570  RAT( x, VAT2(iz, 1, lev)),
571  RAT(tru, VAT2(iz, 1, lev)));
572  } else if (*istop == 3) {
573  Vxcopy(&nxf, &nyf, &nzf, RAT(tru, VAT2(iz, 1, lev)), w1);
574  alpha = -1.0;
575  Vxaxpy(&nxf, &nyf, &nzf, &alpha, RAT(x, VAT2(iz, 1, lev)), w1);
576  rsnrm = Vxnrm2(&nxf, &nyf, &nzf,w1);
577  } else if (*istop == 4) {
578  Vxcopy(&nxf, &nyf, &nzf, RAT(tru, VAT2(iz, 1, lev)), w1);
579  alpha = -1.0;
580  Vxaxpy(&nxf, &nyf, &nzf, &alpha, RAT(x, VAT2(iz, 1, lev)), w1);
581  rsnrm = Vxnrm2(&nxf, &nyf, &nzf,w1);
582  } else if (*istop == 5) {
583  Vxcopy(&nxf, &nyf, &nzf, RAT(tru, VAT2(iz, 1, lev)), w1);
584  alpha = -1.0;
585  Vxaxpy(&nxf, &nyf, &nzf, &alpha, RAT(x, VAT2(iz, 1, lev)), w1);
586  Vnmatvec(&nxf, &nyf, &nzf,
587  RAT(ipc, VAT2(iz, 5, lev)), RAT(rpc, VAT2(iz, 6, lev)),
588  RAT( ac, VAT2(iz, 7, lev)), RAT( cc, VAT2(iz, 1, lev)),
589  w1, w2, w3);
590  rsnrm = VSQRT(Vxdot(&nxf, &nyf, &nzf,w1,w2));
591  } else {
592  VABORT_MSG1("Bad istop value: %d", *istop);
593  }
594 
595  Vprtstp(*iok, *iters, rsnrm, rsden, orsnrm);
596 
597  if ((rsnrm / rsden) <= *errtol)
598  break;
599  }
600 
601  if (*iters >= *itmax) {
602  *ierror = 1;
603  break;
604  }
605  }
606 }
VPUBLIC void Vfmvfas(int *nx, int *ny, int *nz, double *x, int *iz, double *w0, double *w1, double *w2, double *w3, double *w4, int *istop, int *itmax, int *iters, int *ierror, int *nlev, int *ilev, int *nlev_real, int *mgsolv, int *iok, int *iinfo, double *epsiln, double *errtol, double *omega, int *nu1, int *nu2, int *mgsmoo, int *ipc, double *rpc, double *pc, double *ac, double *cc, double *fc, double *tru)
Multigrid nonlinear solve iteration routine.
Definition: mgfasd.c:57
VPUBLIC double Vxdot(int *nx, int *ny, int *nz, double *x, double *y)
Inner product operation for a grid function with boundary values.
Definition: mikpckd.c:173
VEXTERNC void Vnsmooth(int *nx, int *ny, int *nz, int *ipc, double *rpc, double *ac, double *cc, double *fc, double *x, double *w1, double *w2, double *r, int *itmax, int *iters, double *errtol, double *omega, int *iresid, int *iadjoint, int *meth)
call the appropriate non-linear smoothing routine.
Definition: smoothd.c:98
VPUBLIC void Vmvfas(int *nx, int *ny, int *nz, double *x, int *iz, double *w0, double *w1, double *w2, double *w3, double *w4, int *istop, int *itmax, int *iters, int *ierror, int *nlev, int *ilev, int *nlev_real, int *mgsolv, int *iok, int *iinfo, double *epsiln, double *errtol, double *omega, int *nu1, int *nu2, int *mgsmoo, int *ipc, double *rpc, double *pc, double *ac, double *cc, double *fc, double *tru)
Nonlinear multilevel method.
Definition: mgfasd.c:157
VPUBLIC void Vxcopy(int *nx, int *ny, int *nz, double *x, double *y)
A collection of useful low-level routines (timing, etc).
Definition: mikpckd.c:57
VPUBLIC void VinterpPMG(int *nxc, int *nyc, int *nzc, int *nxf, int *nyf, int *nzf, double *xin, double *xout, double *pc)
Apply the prolongation operator.
Definition: matvecd.c:915
VPUBLIC void Vxaxpy(int *nx, int *ny, int *nz, double *alpha, double *x, double *y)
saxpy operation for a grid function with boundary values.
Definition: mikpckd.c:112
VPUBLIC double Vxnrm2(int *nx, int *ny, int *nz, double *x)
Norm operation for a grid function with boundary values.
Definition: mikpckd.c:152
VPUBLIC void Vnmresid(int *nx, int *ny, int *nz, int *ipc, double *rpc, double *ac, double *cc, double *fc, double *x, double *r, double *w1)
Break the matrix data-structure into diagonals and then call the residual routine.
Definition: matvecd.c:598
VPUBLIC void Vazeros(int *nx, int *ny, int *nz, double *x)
Zero out operation for a grid function, including boundary values.
Definition: mikpckd.c:195
VPUBLIC double Vxnrm1(int *nx, int *ny, int *nz, double *x)
Norm operation for a grid function with boundary values.
Definition: mikpckd.c:131
VPUBLIC void Vrestrc(int *nxf, int *nyf, int *nzf, int *nxc, int *nyc, int *nzc, double *xin, double *xout, double *pc)
Apply the restriction operator.
Definition: matvecd.c:782
VPUBLIC void Vextrac(int *nxf, int *nyf, int *nzf, int *nxc, int *nyc, int *nzc, double *xin, double *xout)
Simple injection of a fine grid function into coarse grid.
Definition: matvecd.c:1078
VEXTERNC void Vnmatvec(int *nx, int *ny, int *nz, int *ipc, double *rpc, double *ac, double *cc, double *x, double *y, double *w1)
Break the matrix data-structure into diagonals and then call the matrix-vector routine.
Definition: matvecd.c:232