52 VPUBLIC
void Vmgdriv(
int* iparm,
double* rparm,
53 int* iwork,
double* rwork,
double* u,
54 double* xf,
double* yf,
double* zf,
55 double* gxcf,
double* gycf,
double* gzcf,
56 double* a1cf,
double* a2cf,
double* a3cf,
57 double* ccf,
double* fcf,
double* tcf) {
103 nrwk = VAT(iparm, 1);
104 niwk = VAT(iparm, 2);
108 nlev = VAT(iparm, 6);
111 VASSERT_MSG1(
nlev > 0,
"nlev must be positive: %d",
nlev);
112 VASSERT_MSG1(
nx > 0,
"nx must be positive: %d",
nx);
113 VASSERT_MSG1(
ny > 0,
"nv must be positive: %d",
ny);
114 VASSERT_MSG1(
nz > 0,
"nz must be positive: %d",
nz);
116 mxlv = Vmaxlev(
nx,
ny,
nz);
119 "number of levels exceeds maximum: %d > %d",
141 "real workspace exceeds maximum size: %d > %d",
146 "integer workspace exceeds maximum size: %d > %d",
156 k_cc = k_rpc +
n_rpc;
159 k_ac = k_pc + 27 *
narrc;
164 iz = RAT(iwork, k_iz);
165 ipc = RAT(iwork, k_ipc);
167 rpc = RAT(rwork, k_rpc);
168 pc = RAT(rwork, k_pc);
169 ac = RAT(rwork, k_ac);
170 cc = RAT(rwork, k_cc);
171 fc = RAT(rwork, k_fc);
186 int *
nx,
int *
ny,
int *
nz,
188 int *iz,
int *ipc,
double *rpc,
189 double *pc,
double *ac,
double *cc,
double *fc,
190 double *xf,
double *yf,
double *zf,
191 double *gxcf,
double *gycf,
double *gzcf,
192 double *a1cf,
double *a2cf,
double *a3cf,
193 double *ccf,
double *fcf,
double *tcf) {
227 double tsetupf = 0.0;
228 double tsetupc = 0.0;
239 double errtol_p = 0.0;
241 double rho_min = 0.0;
242 double rho_max = 0.0;
243 double rho_min_mod = 0.0;
244 double rho_max_mod = 0.0;
261 int nlev = VAT(iparm, 6);
267 mgkey = VAT(iparm, 9);
268 itmax = VAT(iparm, 10);
269 istop = VAT(iparm, 11);
270 iinfo = VAT(iparm, 12);
271 ipkey = VAT(iparm, 14);
272 mode = VAT(iparm, 16);
278 iperf = VAT(iparm, 22);
286 Vprtstp(0, -99, 0.0, 0.0, 0.0);
292 Vnm_tstart(30,
"Vmgdrv2: fine problem setup");
299 ipc, rpc, pc, ac, cc, fc,
306 Vnm_tstop(30,
"Vmgdrv2: fine problem setup");
309 Vnm_tstart(30,
"Vmgdrv2: coarse problem setup");
316 ipc, rpc, pc, ac, cc, fc,
323 Vnm_tstop(30,
"Vmgdrv2: coarse problem setup");
326 epsiln = Vnm_epsmac();
345 for (level=1; level <= nlev_real; level++) {
346 nlevd = nlev_real - level + 1;
353 Vmkcors(&numlev, &nxf, &nyf, &nzf, &
nxc, &
nyc, &
nzc);
362 VMESSAGE3(
"Analysis ==> (%3d, %3d, %3d)", nxf, nyf, nzf);
370 VMESSAGE0(
"Power calculating rho(A)");
380 a1cf, a2cf, a3cf, ccf,
381 &rho_max, &rho_max_mod, &errtol_p,
382 &itmax_p, &iters_p, &iinfo_p);
385 VMESSAGE1(
"Power iters = %d", iters_p);
386 VMESSAGE1(
"Power eigmax = %f", rho_max);
387 VMESSAGE1(
"Power (MODEL) = %f", rho_max_mod);
392 VMESSAGE0(
"Ipower calculating lambda_min(A)...");
401 Vipower(&nxf, &nyf, &nzf, u, iz,
402 a1cf, a2cf, a3cf, ccf, fcf,
403 &rho_min, &rho_min_mod, &errtol_p, &itmax_p, &iters_p,
404 &nlevd, &level, &nlev_real, &
mgsolv,
407 ipc, rpc, pc, ac, cc, tcf);
410 VMESSAGE1(
"Ipower iters = %d", iters_p);
411 VMESSAGE1(
"Ipower eigmin = %f", rho_min);
412 VMESSAGE1(
"Ipower (MODEL) = %f", rho_min_mod);
415 VMESSAGE1(
"Condition number = %f", rho_max / rho_min);
416 VMESSAGE1(
"Condition (MODEL) = %f", rho_max_mod / rho_min_mod);
425 VMESSAGE0(
"Mpower calculating rho(M)");
432 Vazeros(&nxf, &nyf, &nzf, RAT(u, VAT2(iz, 1, level)));
435 Vmpower(&nxf, &nyf, &nzf, u, iz,
436 a1cf, a2cf, a3cf, ccf, fcf,
437 &rho_p, &errtol_p, &itmax_p, &iters_p,
438 &nlevd, &level, &nlev_real, &
mgsolv,
439 &iok_p, &iinfo_p, &epsiln,
441 ipc, rpc, pc, ac, cc, fc, tcf);
444 VMESSAGE1(
"Mpower iters = %d", iters_p);
445 VMESSAGE1(
"Mpower rho(M) = %f", rho_p);
451 Vazeros(&nxf, &nyf, &nzf, RAT(u, VAT2(iz, 1, level)));
466 VMESSAGE0(
"Generating algebraic RHS from your soln...");
472 ipc, rpc, ac, cc, ccf, tcf, fc, fcf);
481 Vnm_tstart(30,
"Vmgdrv2: solve");
484 if (mode == 0 || mode == 2) {
492 u, iz, a1cf, a2cf, a3cf, ccf,
494 &ilev, &nlev_real, &
mgsolv,
497 ipc, rpc, pc, ac, cc, fc, tcf);
499 }
else if (
mgkey == 1) {
502 u, iz, a1cf, a2cf, a3cf, ccf,
504 &ilev, &nlev_real, &
mgsolv,
507 ipc, rpc, pc, ac, cc, fc, tcf);
510 VABORT_MSG1(
"Bad mgkey given: %d",
mgkey);
514 if (mode == 1 || mode == 2) {
523 u, iz, a1cf, a2cf, a3cf, ccf, fcf,
525 &ilev, &nlev_real, &
mgsolv,
528 ipc, rpc, pc, ac, cc, fc, tcf);
530 }
else if (
mgkey == 1) {
534 a1cf, a2cf, a3cf, ccf, fcf,
536 &ilev, &nlev_real, &
mgsolv,
539 ipc, rpc, pc, ac, cc, fc, tcf);
542 VABORT_MSG1(
"Bad mgkey given: %d",
mgkey);
547 Vnm_tstop(30,
"Vmgdrv2: solve");
558 int *
nx,
int *
ny,
int *
nz,
564 int *iretot,
int *iintot) {
572 int nc_band, num_band, n_band;
575 int num_nf_oper, num_narrc_oper;
593 for (level=2; level<=*
nlev; level++) {
598 Vmkcors(&numlev, &nxf, &nyf, &nzf,
nxc,
nyc,
nzc);
606 *
narr += nxf * nyf * nzf;
614 }
else if (*
mgdisc == 1) {
617 Vnm_print(2,
"Vmgsz: invalid mgdisc parameter: %d\n", *
mgdisc);
623 }
else if (*
mgcoar == 2) {
626 Vnm_print(2,
"Vmgsz: invalid mgcoar parameter: %d\n", *
mgcoar);
632 }
else if (*
mgsolv == 1) {
634 num_band = 1 + (*
nxc - 2) * (*
nyc - 2);
636 num_band = 1 + (*
nxc - 2) * (*
nyc - 2) + (*
nxc - 2) + 1;
638 nc_band = (*
nxc - 2) * (*
nyc - 2) * (*
nzc - 2);
639 n_band = nc_band * num_band;
641 Vnm_print(2,
"Vmgsz: invalid mgsolv parameter: %d\n", *
mgsolv);
648 *iretot = num_narr * *
narr 649 + (num_nf + num_nf_oper) * *
nf 650 + (num_narrc + num_narrc_oper) * *
narrc VPUBLIC void Vmgdriv(int *iparm, double *rparm, int *iwork, double *rwork, double *u, double *xf, double *yf, double *zf, double *gxcf, double *gycf, double *gzcf, double *a1cf, double *a2cf, double *a3cf, double *ccf, double *fcf, double *tcf)
Multilevel solver driver.
VPUBLIC void VfboundPMG00(int *nx, int *ny, int *nz, double *x)
Initialize a grid function to have a zero boundary value.
VPUBLIC void Vmgdriv2(int *iparm, double *rparm, int *nx, int *ny, int *nz, double *u, int *iz, int *ipc, double *rpc, double *pc, double *ac, double *cc, double *fc, double *xf, double *yf, double *zf, double *gxcf, double *gycf, double *gzcf, double *a1cf, double *a2cf, double *a3cf, double *ccf, double *fcf, double *tcf)
Solves the pde using the multi-grid method.
VEXTERNC void Vmvcs(int *nx, int *ny, int *nz, double *x, int *iz, double *w0, double *w1, double *w2, double *w3, 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)
MG helper functions.
VPUBLIC void Vmgsz(int *mgcoar, int *mgdisc, int *mgsolv, int *nx, int *ny, int *nz, int *nlev, int *nxc, int *nyc, int *nzc, int *nf, int *nc, int *narr, int *narrc, int *n_rpc, int *n_iz, int *n_ipc, int *iretot, int *iintot)
This routine computes the required sizes of the real and integer work arrays for the multigrid code....
VPUBLIC void VfboundPMG(int *ibound, int *nx, int *ny, int *nz, double *x, double *gxc, double *gyc, double *gzc)
Initialize a grid function to have a certain boundary value,.
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.
VPUBLIC void Vazeros(int *nx, int *ny, int *nz, double *x)
Zero out operation for a grid function, including boundary values.
VPUBLIC void Vpower(int *nx, int *ny, int *nz, int *iz, int *ilev, int *ipc, double *rpc, double *ac, double *cc, double *w1, double *w2, double *w3, double *w4, double *eigmax, double *eigmax_model, double *tol, int *itmax, int *iters, int *iinfo)
Power methods for eigenvalue estimation.
VPUBLIC void Vipower(int *nx, int *ny, int *nz, double *u, int *iz, double *w0, double *w1, double *w2, double *w3, double *w4, double *eigmin, double *eigmin_model, double *tol, int *itmax, int *iters, 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 *tru)
Standard inverse power method for minimum eigenvalue estimation.
VPUBLIC void Vbuildops(int *nx, int *ny, int *nz, int *nlev, int *ipkey, int *iinfo, int *ido, int *iz, int *mgprol, int *mgcoar, int *mgsolv, int *mgdisc, int *ipc, double *rpc, double *pc, double *ac, double *cc, double *fc, double *xf, double *yf, double *zf, double *gxcf, double *gycf, double *gzcf, double *a1cf, double *a2cf, double *a3cf, double *ccf, double *fcf, double *tcf)
Build operators, boundary arrays, modify affine vectors ido==0: do only fine level ido==1: do only co...
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.
VPUBLIC void Vbuildstr(int *nx, int *ny, int *nz, int *nlev, int *iz)
Build the nexted operator framework in the array iz.