Skip to content

Commit

Permalink
fixed memory leak with FIELD_FREE
Browse files Browse the repository at this point in the history 
Fixed important memory leak that happened when FIELD_FREE was selected.
Also did a general cleanup of used code.
  • Loading branch information
@tiagopereira
tiagopereira committed Jun 22, 2014
1 parent eb7a5cd commit 00c13b8
Show file tree
Hide file tree
Showing 6 changed files with 28 additions and 244 deletions.
62 changes: 0 additions & 62 deletions rh15d_mpi/background_p.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,8 +2,6 @@
Version: rh2.0, 1.5-D plane-parallel
Author: Tiago Pereira (tiago.pereira@nasa.gov)
Last modified: Tue Jan 24 18:14:00 2012 --
-------------------------- ----------RH-- */

/* Driving subroutine for background opacity sources.
Expand Down Expand Up @@ -200,26 +198,9 @@ void init_Background(bool_t init_BRS)
} else
atmos.backgrrecno = (long *) malloc(spectrum.Nspect * sizeof(long));


/* --- Open background file -- ------------- */

/* get file name, with the MPI rank */
/*
sprintf(file_background,"%s_p%d%s", input.background_File , mpi.rank, fext);
if ((atmos.fd_background =
open(file_background, O_RDWR | O_CREAT | O_TRUNC, PERMISSIONS)) == -1) {
sprintf(messageStr, "Unable to open output file %s",
file_background);
Error(ERROR_LEVEL_2, routineName, messageStr);
}
*/


/* --- Initialise netCDF BRS file -- ------------- */
if (bgdat.write_BRS) init_ncdf_BRS();


/* --- Read background files from Kurucz data file -- ------------- */
atmos.Nrlk = 0;
readKuruczLines(input.KuruczData);
Expand Down Expand Up @@ -724,51 +705,8 @@ void Background_p(bool_t write_analyze_output, bool_t equilibria_only)

if (write_analyze_output) {
/* --- Write background record structure -- ------------ */

if (bgdat.write_BRS) writeBRS_ncdf();

/* --- Write out the metals and molecules -- ------------ */

/* Commented out, not a priority now
writeMetals("metals.out"); // Tiago: must be replaced by NetCDF version
writeMolecules(MOLECULAR_CONCENTRATION_FILE);
*/
}
/* --- Clean up but keep H, H2, and active atom and/or molecule
if appropriate, only for last task -- ------------ */

/* Tiago: commenting this big block for now. Gives problems with rh15d_ray,
and not really necessary (only for tidying up)
if (mpi.task == mpi.Ntasks - 1) {
if (atmos.Natom > 1) {
for (n = 1; n < atmos.Natom; n++)
if (!atmos.atoms[n].active && !atmos.hydrostatic &&
input.solve_ne < ITERATION)
freeAtom(&atmos.atoms[n]);
}
if (atmos.Nmolecule > 1) {
for (n = 1; n < atmos.Nmolecule; n++)
if (!atmos.molecules[n].active &&
!atmos.hydrostatic &&
input.solve_ne < ITERATION)) freeMolecule(&atmos.molecules[n]);
}
if (strcmp(input.KuruczData, "none")) {
free(atmos.Tpf); atmos.Tpf = NULL;
for (n = 0; n < atmos.Nelem; n++) {
free(atmos.elements[n].ionpot);
freeMatrix((void **) atmos.elements[n].pf);
if (atmos.elements[n].n)
freeMatrix((void **) atmos.elements[n].n);
}
}
if (do_fudge) {
free(lambda_fudge);
freeMatrix((void **) fudge);
}
} */


getCPU(3, TIME_POLL, "Background Opacity");

Expand Down
54 changes: 3 additions & 51 deletions rh15d_mpi/distribute_jobs.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -42,18 +42,10 @@ extern NCDF_Atmos_file infile;


/* ------- begin -------------------------- distribute_jobs.c --- */
void distribute_jobs(void)
void distribute_jobs(void)
/* Distributes the work between the available processes */
{
/* The end product of this routine must be:
For each process:
* mpi.Ntasks --DONE
* An array of rank (mpi.Ntasks,2) that gives (mpi.ix, mpi.iy)
for each task index --DONE (mpi.my_tasks)
* A map of (mpi.ix, mpi.iy) -> (xi, yi) --DONE (two 1D arrays: mpi.xnum, mpi.ynum)
*/


long *tasks, remain_tasks, i, j;

mpi.backgrrecno = 0;
Expand Down Expand Up @@ -97,50 +89,10 @@ void distribute_jobs(void)

mpi.total_tasks = remain_tasks;

/* Abort if more processes than tasks (avoid idle processes waiting forever)
if (mpi.size > remain_tasks) {
sprintf(messageStr,
"\n*** More MPI processes (%d) than tasks (%ld), aborting.\n",
mpi.size, remain_tasks);
Error(ERROR_LEVEL_2, "distribute_jobs", messageStr);
}
/* Calculate tasks and distribute */
tasks = get_tasks(remain_tasks, mpi.size);
mpi.Ntasks = tasks[mpi.rank];
mpi.taskmap = get_taskmap(remain_tasks, tasks, &mpi.my_start);

// TIAGO TEMPORARY for specific x, y points:
/*
mpi.ynum[0] = 39;
mpi.ynum[1] = 78;
mpi.ynum[2] = 184;
mpi.ynum[3] = 189;
mpi.ynum[4] = 228;
mpi.ynum[5] = 244;
mpi.ynum[6] = 280;
mpi.ynum[7] = 284;
mpi.ynum[8] = 342;
mpi.ynum[9] = 347;
mpi.ynum[10] = 452;
mpi.ynum[11] = 466;
mpi.xnum[0] = 36;
mpi.xnum[1] = 70;
mpi.xnum[2] = 112;
mpi.xnum[3] = 162;
mpi.xnum[4] = 163;
mpi.xnum[5] = 210;
mpi.xnum[6] = 289;
mpi.xnum[7] = 246;
mpi.xnum[8] = 303;
mpi.xnum[9] = 370;
mpi.xnum[10] = 387;
mpi.xnum[11] = 388;
mpi.xnum[12] = 449;
*/

free(tasks);

return;
Expand Down
49 changes: 0 additions & 49 deletions rh15d_mpi/ncdfatmos.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,8 +2,6 @@
Version: rh2.0, 1.5-D plane-parallel
Author: Tiago Pereira (tiago.pereira@nasa.gov)
Last modified: Tue Nov 23 16:03:23 2010 --
-------------------------- ----------RH-- */

/* --- Reads atmospheric model in NetCDF format. -- -------------- */
Expand Down Expand Up @@ -229,8 +227,6 @@ void readAtmos_ncdf(int xi, int yi, Atmosphere *atmos, Geometry *geometry,
} else {
mpi.zcut = 0;
}

//printf("Process %d: zcut = %d\n", mpi.rank, mpi.zcut);

/* Get z again */
start[0] = input.p15d_nt; count[0] = 1;
Expand Down Expand Up @@ -258,35 +254,8 @@ void readAtmos_ncdf(int xi, int yi, Atmosphere *atmos, Geometry *geometry,
if (infile->vturb_varid != -1) {
if ((ierror = nc_get_vara_double(ncid, infile->vturb_varid, &start[3], &count[3],
atmos->vturb))) ERR(ierror,routineName);
} else {
/* Tiago's microturbulence fudge */
/* Find temperature minimum
imin = 0;
i50k = 0;
Tmin = 1000000.;
diff = 1000000.;
for (j = 0; j < atmos->Nspace; j++) {
if (atmos->T[j] < Tmin) {
Tmin = atmos->T[j];
imin = j;
}
if (fabs(atmos->T[j] - 50000.) < diff) {
diff = fabs(atmos->T[j] - 50000.);
i50k = j;
}
}
//printf("imin, i50k = %i, %i T=%e, %e\n", imin, i50k, atmos->T[imin], atmos->T[i50k]);
/* Scale vturb
for (j = 0; j < imin+1; j++) {
atmos->vturb[j] = 1.e3 + 1.e4*pow(atmos->T[j] - Tmin, 0.333)/pow(50000. - Tmin, 0.333);
//printf("z, Temp, vturb = %.3e %.3e %.3e\n", geometry->height[j]/1.e6, atmos->T[j], atmos->vturb[j]/1.e3);
}
*/
}


/* Read magnetic field */
if (atmos->Stokes) {
Bx = (double *) malloc(atmos->Nspace * sizeof(double));
Expand Down Expand Up @@ -542,13 +511,6 @@ void depth_refine(Atmosphere *atmos, Geometry *geometry, double Tmax) {
for (k = 1; k <= k1; k++)
aind[k] *= (atmos->Nspace-1)/aind[k1];

/*
printf("Original quantities:\n---------------------\n");
for (k = 0; k < atmos->Nspace; k++)
printf("%4li %12.4e %12.4e %12.4e %12.4e %12.4e %12.4e\n", k, geometry->height[k],atmos->T[k], atmos->ne[k], geometry->vel[k], atmos->B[k], atmos->gamma_B[k]);
*/


/* --- Create new height scale --- */
splineCoef(k1-k0+1, &aind[k0], &geometry->height[k0]);
splineEval(atmos->Nspace, xpt, new_height, hunt=FALSE);
Expand All @@ -561,10 +523,6 @@ void depth_refine(Atmosphere *atmos, Geometry *geometry, double Tmax) {
atmos->nH[i][k] = log(atmos->nH[i][k]);
}

/*
Linear(atmos->Nspace, geometry->height, atmos->T, atmos->Nspace, new_height,
buf, FALSE);
*/
splineCoef(atmos->Nspace, geometry->height, atmos->T);
splineEval(atmos->Nspace, new_height, buf, hunt=FALSE);
memcpy((void *) atmos->T, (void *) buf, bufsize);
Expand Down Expand Up @@ -610,13 +568,6 @@ void depth_refine(Atmosphere *atmos, Geometry *geometry, double Tmax) {
}

memcpy((void *) geometry->height, (void *) new_height, bufsize);

/*
printf("Interpolated quantities:\n---------------------\n");
for (k = 0; k < atmos->Nspace; k++)
printf("%4li %12.4e %12.4e %12.4e %12.4e %12.4e %12.4e\n", k, geometry->height[k],atmos->T[k], atmos->ne[k], geometry->vel[k],atmos->vturb[k] , 0.);//atmos->B[k], atmos->gamma_B[k]);
*/


if (nhm_flag) {
free(atmos->nHmin);
Expand Down
30 changes: 6 additions & 24 deletions rh15d_mpi/parallel.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,8 +2,6 @@
Version: rh2.0, 1.5-D plane-parallel
Author: Tiago Pereira (tiago.pereira@nasa.gov)
Last modified: Mon Jan 03 14:28:25 2011 --
-------------------------- -----------RH-- */

/* --- Set of tools to deal with IO initialisation and closure ------- */
Expand Down Expand Up @@ -125,7 +123,6 @@ void initParallelIO(bool_t run_ray, bool_t writej) {
/* ------- end -------------------------- initParallelIO.c --- */



/* ------- begin -------------------------- closeParallelIO.c --- */
void closeParallelIO(bool_t run_ray, bool_t writej) {

Expand Down Expand Up @@ -223,8 +220,6 @@ void UpdateAtmosDep(void) {

/* Reallocate some stuff (because of varying Nspace) */



/* Free collision rate array, will be reallocated by calls in Background_p */
if (atom->C != NULL) {
freeMatrix((void **) atom->C);
Expand Down Expand Up @@ -261,7 +256,7 @@ void UpdateAtmosDep(void) {
line->wphi = NULL;
}

if (atmos.moving && line->polarizable && (input.StokesMode>FIELD_FREE)) {
if (atmos.moving && line->polarizable && (input.StokesMode >= FIELD_FREE)) {

if (line->phi_Q != NULL) {
freeMatrix((void **) line->phi_Q);
Expand Down Expand Up @@ -320,8 +315,8 @@ void UpdateAtmosDep(void) {
else
Nlamu = line->Nlambda;

// Idea: instead of doing this, why not free and just set line->rho_prd = NULL,
// (and also line->Qelast?), because profile.c will act on that and reallocate
/* Idea: instead of doing this, why not free and just set line->rho_prd = NULL,
(and also line->Qelast?), because profile.c will act on that and reallocate */
if (line->rho_prd != NULL) freeMatrix((void **) line->rho_prd);
line->rho_prd = matrix_double(Nlamu, atmos.Nspace);

Expand Down Expand Up @@ -507,12 +502,6 @@ void copyBufVars(bool_t writej) {
AtomicLine *line;
AtomicContinuum *continuum;


/* J, J20
memcpy((void *) &iobuf.J[ind*spectrum.Nspect], (void *) spectrum.J[0],
spectrum.Nspect * atmos.Nspace * sizeof(double));
*/

if (writej) {
for (nspect=0; nspect < spectrum.Nspect; nspect++) {
i = 0;
Expand Down Expand Up @@ -571,10 +560,8 @@ void copyBufVars(bool_t writej) {
molecule->Nv * atmos.Nspace * sizeof(double));

}

ind += atmos.Nspace;


return;
}

Expand Down Expand Up @@ -723,9 +710,9 @@ void writeOutput(bool_t writej) {
MPI_Status status;

/* Write output in order of rank. First 0, then send to next, until all
processes have written the output. */
processes have written the output.
/* This can also be used with an integer, like if mpi.rank > ml,
This can also be used with an integer, like if mpi.rank > ml,
and then adding ml to mpi.rank in the send. But for now not using
if (mpi.rank > 0)
MPI_Recv(&msg, 0, MPI_INT, MPI_ANY_SOURCE, MPI_ANY_TAG, mpi.comm, &status);
Expand All @@ -743,17 +730,12 @@ void writeOutput(bool_t writej) {
writeMPI_all();
if (writej) writeJ_all();
writeAux_all();
//writeAtmos_all();
/* writeAtmos_all(); */

sprintf(messageStr, "Process %3d: *** END output\n", mpi.rank);
fprintf(mpi.main_logfile, messageStr);
Error(MESSAGE, "main", messageStr);
}
/*
if (mpi.rank < mpi.size - 1) {
MPI_Send(0, 0, MPI_INT, mpi.rank + 1, 111, mpi.comm);
}
*/



Expand Down
13 changes: 2 additions & 11 deletions rh15d_mpi/rh15d_ray.c
View file
Original file line number Diff line number Diff line change
@@ -1,4 +1,3 @@
/* $Id$ */
#include <string.h>
#include <stdlib.h>
#include <math.h>
Expand Down Expand Up @@ -53,18 +52,11 @@ IO_buffer iobuf;
int main(int argc, char *argv[])
{
bool_t write_analyze_output, equilibria_only, run_ray, writej, exit_on_EOF;
int niter, nact, i, Nspect, Ntest, k, Nread, Nrequired, ierror,
checkPoint, save_Nrays, *wave_index = NULL, conv_iter;
int niter, i, Nspect, Nread, Nrequired,
checkPoint, save_Nrays, *wave_index = NULL;
double muz, save_muz, save_mux, save_muy, save_wmu;
Atom *atom;
Molecule *molecule;
AtomicLine *line;
ActiveSet *as;
FILE *fp_ray;

time_t curtime;
struct tm *loctime;
char timestr[ARR_STRLEN];
char inputLine[MAX_LINE_SIZE];

/* --- Set up MPI ---------------------- -------------- */
Expand Down Expand Up @@ -249,7 +241,6 @@ int main(int argc, char *argv[])

calculate_ray();
writeRay();


/* Put back previous values for geometry */
atmos.Nrays = geometry.Nrays = save_Nrays;
Expand Down
Loading

0 comments on commit 00c13b8

Please sign in to comment.