APBS 3.0.0
Loading...
Searching...
No Matches
vpmgp.c
Go to the documentation of this file.
1
57#include "vpmgp.h"
58
59VEMBED(rcsid="$Id$")
60
61/* ///////////////////////////////////////////////////////////////////////////
62// Class Vpmgp: Inlineable methods
64#if !defined(VINLINE_VACC)
65#endif /* if !defined(VINLINE_VACC) */
66
67/* ///////////////////////////////////////////////////////////////////////////
68// Class Vpmgp: Non-inlineable methods
70
71/* ///////////////////////////////////////////////////////////////////////////
72// Routine: Vpmgp_ctor
73//
74// Author: Nathan Baker
76VPUBLIC Vpmgp* Vpmgp_ctor(MGparm *mgparm) {
77
78 Vpmgp *thee = VNULL;
79
80 /* Set up the structure */
81 thee = (Vpmgp*)Vmem_malloc(VNULL, 1, sizeof(Vpmgp) );
82 VASSERT( thee != VNULL);
83 VASSERT(Vpmgp_ctor2(thee,mgparm));
84
85 return thee;
86}
87
88/* ///////////////////////////////////////////////////////////////////////////
89// Routine: Vpmgp_ctor2
90//
91// Author: Nathan Baker
93VPUBLIC int Vpmgp_ctor2(Vpmgp *thee,MGparm *mgparm) {
94
95 /* Specified parameters */
96 thee->nx = mgparm->dime[0];
97 thee->ny = mgparm->dime[1];
98 thee->nz = mgparm->dime[2];
99 thee->hx = mgparm->grid[0];
100 thee->hy = mgparm->grid[1];
101 thee->hzed = mgparm->grid[2];
102 thee->xlen = ((double)(mgparm->dime[0]-1))*mgparm->grid[0];
103 thee->ylen = ((double)(mgparm->dime[1]-1))*mgparm->grid[1];
104 thee->zlen = ((double)(mgparm->dime[2]-1))*mgparm->grid[2];
105 thee->nlev = mgparm->nlev;
106
107 thee->nonlin = mgparm->nonlintype;
108 thee->meth = mgparm->method;
109
110#ifdef DEBUG_MAC_OSX_OCL
111#include "mach_chud.h"
112 if(kOpenCLAvailable)
113 thee->meth = 4;
114#endif
115
116 if (thee->nonlin == NONLIN_LPBE) thee->ipkey = IPKEY_LPBE; /* LPBE case */
117 else if(thee->nonlin == NONLIN_SMPBE) thee->ipkey = IPKEY_SMPBE; /* SMPBE case */
118 else thee->ipkey = IPKEY_NPBE; /* NPBE standard case */
119
120 /* Default parameters */
121 if (mgparm->setetol) { /* If etol is set by the user in APBS input file, then use this custom-defined etol */
122 thee->errtol = mgparm->etol;
123 Vnm_print(1, " Error tolerance (etol) is now set to user-defined \
124value: %g \n", thee->errtol);
125 Vnm_print(0, "Error tolerance (etol) is now set to user-defined \
126value: %g \n", thee->errtol);
127 } else thee->errtol = 1.0e-6; /* Here are a few comments. Mike had this set to
128 * 1e-9; convential wisdom sets this at 1e-6 for
129 * the PBE; Ray Luo sets this at 1e-3 for his
130 * accelerated PBE (for dynamics, etc.) */
131 thee->itmax = 200;
132 thee->istop = 1;
133 thee->iinfo = 1; /* I'd recommend either 1 (for debugging LPBE) or 2 (for debugging NPBE), higher values give too much output */
134
135 thee->bcfl = BCFL_SDH;
136 thee->key = 0;
137 thee->iperf = 0;
138 thee->mgcoar = 2;
139 thee->mgkey = 0;
140 thee->nu1 = 2;
141 thee->nu2 = 2;
142 thee->mgprol = 0;
143 thee->mgdisc = 0;
144 thee->omegal = 19.4e-1;
145 thee->omegan = 9.0e-1;
146 thee->ipcon = 3;
147 thee->irite = 8;
148 thee->xcent = 0.0;
149 thee->ycent = 0.0;
150 thee->zcent = 0.0;
151
152 /* Default value for all APBS runs */
153 thee->mgsmoo = 1;
154 if (thee->nonlin == NONLIN_NPBE || thee->nonlin == NONLIN_SMPBE) {
155 /* SMPBE Added - SMPBE needs to mimic NPBE */
156 Vnm_print(0, "Vpmp_ctor2: Using meth = 1, mgsolv = 0\n");
157 thee->mgsolv = 0;
158 } else {
159 /* Most rigorous (good for testing) */
160 Vnm_print(0, "Vpmp_ctor2: Using meth = 2, mgsolv = 1\n");
161 thee->mgsolv = 1;
162 }
163
164 /* TEMPORARY USEAQUA */
165 /* If we are using aqua, our solution method is either VSOL_CGMGAqua or VSOL_NewtonAqua
166 * so we need to temporarily override the mgsolve value and set it to 0
167 */
168 if(mgparm->useAqua == 1) thee->mgsolv = 0;
169
170 return 1;
171}
172
173/* ///////////////////////////////////////////////////////////////////////////
174// Routine: Vpmgp_dtor
175//
176// Author: Nathan Baker
178VPUBLIC void Vpmgp_dtor(Vpmgp **thee) {
179
180 if ((*thee) != VNULL) {
181 Vpmgp_dtor2(*thee);
182 Vmem_free(VNULL, 1, sizeof(Vpmgp), (void **)thee);
183 (*thee) = VNULL;
184 }
185
186}
187
188/* ///////////////////////////////////////////////////////////////////////////
189// Routine: Vpmgp_dtor2
190//
191// Author: Nathan Baker
193VPUBLIC void Vpmgp_dtor2(Vpmgp *thee) { ; }
194
195
196VPUBLIC void Vpmgp_size(
197 Vpmgp *thee
198 )
199{
200
201 int num_nf = 0;
202 int num_narr = 2;
203 int num_narrc = 27;
204 int nxf, nyf, nzf, level, num_nf_oper, num_narrc_oper, n_band, nc_band, num_band, iretot;
205
206 thee->nf = thee->nx * thee->ny * thee->nz;
207 thee->narr = thee->nf;
208 nxf = thee->nx;
209 nyf = thee->ny;
210 nzf = thee->nz;
211 thee->nxc = thee->nx;
212 thee->nyc = thee->ny;
213 thee->nzc = thee->nz;
214
215 for (level=2; level<=thee->nlev; level++) {
216 Vpmgp_makeCoarse(1, nxf, nyf, nzf, &(thee->nxc), &(thee->nyc), &(thee->nzc)); /* NAB TO-DO -- implement this function and check which variables need to be passed by reference... */
217 nxf = thee->nxc;
218 nyf = thee->nyc;
219 nzf = thee->nzc;
220 thee->narr = thee->narr + (nxf * nyf * nzf);
221 }
222
223 thee->nc = thee->nxc * thee->nyc * thee->nzc;
224 thee->narrc = thee->narr - thee->nf;
225
226 /* Box or FEM discretization on fine grid? */
227 switch (thee->mgdisc) { /* NAB TO-DO: This needs to be changed into an enumeration */
228 case 0:
229 num_nf_oper = 4;
230 break;
231 case 1:
232 num_nf_oper = 14;
233 break;
234 default:
235 Vnm_print(2, "Vpmgp_size: Invalid mgdisc value (%d)!\n", thee->mgdisc);
236 VASSERT(0);
237 }
238
239 /* Galerkin or standard coarsening? */
240 switch (thee->mgcoar) { /* NAB TO-DO: This needs to be changed into an enumeration */
241 case 0:
242 if (thee->mgdisc != 0) {
243 Vnm_print(2, "Vpmgp_size: Invalid mgcoar value (%d); must be used with mgdisc 0!\n", thee->mgcoar);
244 VASSERT(0);
245 }
246 num_narrc_oper = 4;
247 break;
248 case 1:
249 if (thee->mgdisc != 0) {
250 Vnm_print(2, "Vpmgp_size: Invalid mgcoar value (%d); must be used with mgdisc 0!\n", thee->mgcoar);
251 VASSERT(0);
252 }
253 num_narrc_oper = 14;
254 break;
255 case 2:
256 num_narrc_oper = 14;
257 break;
258 default:
259 Vnm_print(2, "Vpmgp_size: Invalid mgcoar value (%d)!\n", thee->mgcoar);
260 VASSERT(0);
261 }
262
263 /* LINPACK storage on coarse grid */
264 switch (thee->mgsolv) { /* NAB TO-DO: This needs to be changed into an enumeration */
265 case 0:
266 n_band = 0;
267 break;
268 case 1:
269 if ( ( (thee->mgcoar == 0) || (thee->mgcoar == 1)) && (thee->mgdisc == 0) ) {
270 num_band = 1 + (thee->nxc-2)*(thee->nyc-2);
271 } else {
272 num_band = 1 + (thee->nxc-2)*(thee->nyc-2) + (thee->nxc-2) + 1;
273 }
274 nc_band = (thee->nxc-2)*(thee->nyc-2)*(thee->nzc-2);
275 n_band = nc_band * num_band;
276 break;
277 default:
278 Vnm_print(2, "Vpmgp_size: Invalid mgsolv value (%d)!\n", thee->mgsolv);
279 VASSERT(0);
280 }
281
282 /* Real storage parameters */
283 thee->n_rpc = 100*(thee->nlev+1);
284
285 /* Resulting total required for real storage */
286 thee->nrwk = num_narr*thee->narr + (size_t)(num_nf + num_nf_oper)*thee->nf + (size_t)(num_narrc + num_narrc_oper)*thee->narrc + n_band + thee->n_rpc;
287
288 /* Integer storage parameters */
289 thee->n_iz = 50*(thee->nlev+1);
290 thee->n_ipc = 100*(thee->nlev+1);
291 thee->niwk = thee->n_iz + thee->n_ipc;
292}
293
294VPRIVATE int coarsenThis(int nOld) {
295
296 int nOut;
297
298 nOut = (nOld - 1) / 2 + 1;
299
300 if (((nOut-1)*2) != (nOld-1)) {
301 Vnm_print(2, "Vpmgp_makeCoarse: Warning! The grid dimensions you have chosen are not consistent with the nlev you have specified!\n");
302 Vnm_print(2, "Vpmgp_makeCoarse: This calculation will only work if you are running with mg-dummy type.\n");
303 }
304 if (nOut < 1) {
305 Vnm_print(2, "D'oh! You coarsened the grid below zero! How did you do that?\n");
306 VASSERT(0);
307 }
308
309 return nOut;
310}
311
312VPUBLIC void Vpmgp_makeCoarse(
313 int numLevel,
314 int nxOld,
315 int nyOld,
316 int nzOld,
317 int *nxNew,
318 int *nyNew,
319 int *nzNew
320 )
321{
322 int nxtmp, nytmp, nztmp, iLevel;
323
324 for (iLevel=0; iLevel<numLevel; iLevel++) {
325 nxtmp = *nxNew;
326 nytmp = *nyNew;
327 nztmp = *nzNew;
328 *nxNew = coarsenThis(nxtmp);
329 *nyNew = coarsenThis(nytmp);
330 *nzNew = coarsenThis(nztmp);
331 }
332
333
334}
#define VEMBED(rctag)
Allows embedding of RCS ID tags in object files.
Definition vhal.h:556
@ BCFL_SDH
Definition vhal.h:209
@ IPKEY_LPBE
Definition vhal.h:159
@ IPKEY_NPBE
Definition vhal.h:160
@ IPKEY_SMPBE
Definition vhal.h:158
VPUBLIC void Vpmgp_makeCoarse(int numLevel, int nxOld, int nyOld, int nzOld, int *nxNew, int *nyNew, int *nzNew)
Coarsen the grid by the desired number of levels and determine the resulting numbers of grid points.
Definition vpmgp.c:312
VPUBLIC int Vpmgp_ctor2(Vpmgp *thee, MGparm *mgparm)
FORTRAN stub to construct PMG parameter object and initialize to default values.
Definition vpmgp.c:93
Contains public data members for Vpmgp class/module.
Definition vpmgp.h:80
Contains declarations for class Vpmgp.