Skip to content

Commit

Permalink
Updated testing for better performance when recreating all tests. Als…
Browse files Browse the repository at this point in the history 
…o rebuilt all tests to account for new outputs. Fixed issue where last value of DeltaTime in log did not equal last value in output file.
  • Loading branch information
@RoryBarnes
RoryBarnes committed Feb 5, 2024
1 parent ae3a1e1 commit d15e244
Show file tree
Hide file tree
Showing 53 changed files with 31,824 additions and 16,039 deletions.
324 changes: 185 additions & 139 deletions src/atmesc.c
View file

Large diffs are not rendered by default.

56 changes: 28 additions & 28 deletions src/eqtide.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1490,7 +1490,8 @@ void VerifyCPL(BODY *body, CONTROL *control, FILES *files, OPTIONS *options,
/* Tidal Q was also set */
if (control->Io.iVerbose >= VERBINPUT) {
fprintf(stderr,
"WARNING: Phase lag model selected, but both %s and %s set in file %s. "
"WARNING: Phase lag model selected, but both %s and %s set in "
"file %s. "
"Using %s = %lf and ignoring %s.\n",
options[OPT_TIDALTAU].cName, options[OPT_TIDALQ].cName,
options[OPT_TIDALTAU].cFile[iBody + 1],
Expand Down Expand Up @@ -1591,10 +1592,9 @@ void VerifyCPL(BODY *body, CONTROL *control, FILES *files, OPTIONS *options,
void VerifyPerturbersEqtide(BODY *body, FILES *files, OPTIONS *options,
UPDATE *update, int iNumBodies, int iBody) {
int iPert, iBodyPert, iVar, ok;
int *bFound = malloc(iNumBodies*sizeof(int));
int *bFound = malloc(iNumBodies * sizeof(int));

for (iBody = 0; iBody < iNumBodies; iBody++) {
fprintf(stderr,"Body: %s\n",body[iBody].cName);
if (body[iBody].bEqtide) {
if (body[iBody].iTidePerts > 0) {
for (iPert = 0; iPert < body[iBody].iTidePerts; iPert++) {
Expand Down Expand Up @@ -2665,8 +2665,8 @@ void WriteEqRotPer(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system,
UNITS *units, UPDATE *update, int iBody, double *dTmp,
char cUnit[]) {

int iOrbiter = fiAssignTidalOrbiter(body,iBody);
*dTmp = fdFreqToPer(fdEqRotRate(
int iOrbiter = fiAssignTidalOrbiter(body, iBody);
*dTmp = fdFreqToPer(fdEqRotRate(
body, iBody, body[iOrbiter].dMeanMotion, body[iOrbiter].dEccSq,
control->Evolve.iEqtideModel, control->Evolve.bDiscreteRot));

Expand All @@ -2683,7 +2683,7 @@ void WriteEqRotPerCont(BODY *body, CONTROL *control, OUTPUT *output,
SYSTEM *system, UNITS *units, UPDATE *update, int iBody,
double *dTmp, char cUnit[]) {

int iOrbiter = fiAssignTidalOrbiter(body,iBody);
int iOrbiter = fiAssignTidalOrbiter(body, iBody);

// To CPL, or to CTL? That is the question
if (control->Evolve.iEqtideModel == CPL) {
Expand All @@ -2708,7 +2708,7 @@ void WriteEqRotPerDiscrete(BODY *body, CONTROL *control, OUTPUT *output,
SYSTEM *system, UNITS *units, UPDATE *update,
int iBody, double *dTmp, char cUnit[]) {

int iOrbiter = fiAssignTidalOrbiter(body,iBody);
int iOrbiter = fiAssignTidalOrbiter(body, iBody);
if (control->Evolve.iEqtideModel == CPL) {
*dTmp = fdFreqToPer(fdCPLEqRotRateDiscrete(body[iOrbiter].dMeanMotion,
body[iOrbiter].dEccSq));
Expand All @@ -2729,10 +2729,10 @@ void WriteEqRotRate(BODY *body, CONTROL *control, OUTPUT *output,
SYSTEM *system, UNITS *units, UPDATE *update, int iBody,
double *dTmp, char cUnit[]) {

int iOrbiter = fiAssignTidalOrbiter(body,iBody);
*dTmp = fdEqRotRate(body, iBody, body[iOrbiter].dMeanMotion,
body[iOrbiter].dEccSq, control->Evolve.iEqtideModel,
control->Evolve.bDiscreteRot);
int iOrbiter = fiAssignTidalOrbiter(body, iBody);
*dTmp = fdEqRotRate(body, iBody, body[iOrbiter].dMeanMotion,
body[iOrbiter].dEccSq, control->Evolve.iEqtideModel,
control->Evolve.bDiscreteRot);

if (output->bDoNeg[iBody]) {
*dTmp *= output->dNeg;
Expand All @@ -2747,7 +2747,7 @@ void WriteEqRotRateCont(BODY *body, CONTROL *control, OUTPUT *output,
SYSTEM *system, UNITS *units, UPDATE *update, int iBody,
double *dTmp, char cUnit[]) {

int iOrbiter = fiAssignTidalOrbiter(body,iBody);
int iOrbiter = fiAssignTidalOrbiter(body, iBody);

// To CPL, or to CTL? That is the question
if (control->Evolve.iEqtideModel == CPL) {
Expand All @@ -2770,7 +2770,7 @@ void WriteEqRotRateCont(BODY *body, CONTROL *control, OUTPUT *output,
void WriteEqRotRateDiscrete(BODY *body, CONTROL *control, OUTPUT *output,
SYSTEM *system, UNITS *units, UPDATE *update,
int iBody, double *dTmp, char cUnit[]) {
int iOrbiter = fiAssignTidalOrbiter(body,iBody);
int iOrbiter = fiAssignTidalOrbiter(body, iBody);

if (control->Evolve.iEqtideModel == CPL) {
*dTmp = fdCPLEqRotRateDiscrete(body[iOrbiter].dMeanMotion,
Expand All @@ -2792,7 +2792,7 @@ void WriteEqTidePower(BODY *body, CONTROL *control, OUTPUT *output,
SYSTEM *system, UNITS *units, UPDATE *update, int iBody,
double *dTmp, char cUnit[]) {

int iOrbiter = fiAssignTidalOrbiter(body,iBody);
int iOrbiter = fiAssignTidalOrbiter(body, iBody);
// Why is this Eq??????? XXX
if (control->Evolve.iEqtideModel == CPL) {
*dTmp = fdCPLTidePowerEq(body[iBody].dTidalZ[iOrbiter], body[iBody].dEccSq,
Expand Down Expand Up @@ -2932,13 +2932,13 @@ void WritePowerEqtide(BODY *body, CONTROL *control, OUTPUT *output,
//*dTmp = body[iBody].dTidalPowMan;


int jBody=0;
if (iBody==0)
int jBody = 0;
if (iBody == 0)
jBody++;
/*
fprintf(stderr,"\niBody: %d\n",iBody);
fprintf(stderr,"\njBody: %d\n",jBody);
fprintf(stderr,"TidalZ: %lf\n",body[iBody].dTidalZ[jBody]);
/*
fprintf(stderr,"\niBody: %d\n",iBody);
fprintf(stderr,"\njBody: %d\n",jBody);
fprintf(stderr,"TidalZ: %lf\n",body[iBody].dTidalZ[jBody]);
*/

if (output->bDoNeg[iBody]) {
Expand Down Expand Up @@ -3433,7 +3433,7 @@ void LogBodyEqtide(BODY *body, CONTROL *control, OUTPUT *output, SYSTEM *system,
fprintf(fp, "----- EQTIDE PARAMETERS (%s)------\n", body[iBody].cName);
for (iOut = iStart; iOut < OUTENDEQTIDE; iOut++) {
if (output[iOut].iNum > 0) {
//fprintf(stderr,"iOut = %d.\n",iOut);
// fprintf(stderr,"iOut = %d.\n",iOut);
WriteLogEntry(body, control, &output[iOut], system, update, fnWrite[iOut],
fp, iBody);
}
Expand Down Expand Up @@ -3906,11 +3906,11 @@ void ForceBehaviorEqtide(BODY *body, MODULE *module, EVOLVE *evolve, IO *io,

/**
* Identify and return the index for the perturbing body.
*
*
*/
*
*
*/

int fiAssignTidalPerturber(BODY *body,int iBody) {
int fiAssignTidalPerturber(BODY *body, int iBody) {
int iPerturber;

if (!bPrimary(body, iBody)) {
Expand All @@ -3925,7 +3925,7 @@ int fiAssignTidalPerturber(BODY *body,int iBody) {

/* Non-central body contains the orbital parameters */

int fiAssignTidalOrbiter(BODY *body,int iBody) {
int fiAssignTidalOrbiter(BODY *body, int iBody) {
int iOrbiter;

if (!bPrimary(body, iBody)) {
Expand Down Expand Up @@ -3963,8 +3963,8 @@ double fdCPLTidePower(BODY *body, int iBody) {

// Does this work with DF's changes to da/dt with the synchronous case?
// See Fleming et al., 2018
//fprintf(stderr,"\niBody: %d\n",iBody);
//fprintf(stderr,"TidalZ[%d]: %lf\n",iIndex,body[iBody].dTidalZ[iIndex]);
// fprintf(stderr,"\niBody: %d\n",iBody);
// fprintf(stderr,"TidalZ[%d]: %lf\n",iIndex,body[iBody].dTidalZ[iIndex]);

dOrbPow += -body[iBody].dTidalZ[iIndex] / 8 *
(4 * body[iBody].iTidalEpsilon[iIndex][0] +
Expand Down
34 changes: 21 additions & 13 deletions src/evolve.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -423,7 +423,8 @@ void RungeKutta4Step(BODY *body, CONTROL *control, SYSTEM *system,
iNumEqns = update[iBody].iNumEqns[iVar];
for (iEqn = 0; iEqn < iNumEqns; iEqn++) {
daDerivVar += iDir * evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
evolve->daDerivProc[0][iBody][iVar][iEqn] = evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
evolve->daDerivProc[0][iBody][iVar][iEqn] =
evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
}
evolve->daDeriv[0][iBody][iVar] = daDerivVar;
}
Expand Down Expand Up @@ -465,7 +466,8 @@ void RungeKutta4Step(BODY *body, CONTROL *control, SYSTEM *system,
iNumEqns = update[iBody].iNumEqns[iVar];
for (iEqn = 0; iEqn < iNumEqns; iEqn++) {
daDerivVar += iDir * evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
evolve->daDerivProc[1][iBody][iVar][iEqn] = evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
evolve->daDerivProc[1][iBody][iVar][iEqn] =
evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
}
evolve->daDeriv[1][iBody][iVar] = daDerivVar;
}
Expand Down Expand Up @@ -508,7 +510,8 @@ void RungeKutta4Step(BODY *body, CONTROL *control, SYSTEM *system,
iNumEqns = update[iBody].iNumEqns[iVar];
for (iEqn = 0; iEqn < iNumEqns; iEqn++) {
daDerivVar += iDir * evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
evolve->daDerivProc[2][iBody][iVar][iEqn] = evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
evolve->daDerivProc[2][iBody][iVar][iEqn] =
evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
}
evolve->daDeriv[2][iBody][iVar] = daDerivVar;
}
Expand Down Expand Up @@ -557,7 +560,8 @@ void RungeKutta4Step(BODY *body, CONTROL *control, SYSTEM *system,
iNumEqns = update[iBody].iNumEqns[iVar];
for (iEqn = 0; iEqn < iNumEqns; iEqn++) {
daDerivVar += iDir * evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
evolve->daDerivProc[3][iBody][iVar][iEqn] = evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
evolve->daDerivProc[3][iBody][iVar][iEqn] =
evolve->tmpUpdate[iBody].daDerivProc[iVar][iEqn];
}
evolve->daDeriv[3][iBody][iVar] = daDerivVar;
}
Expand All @@ -574,11 +578,12 @@ void RungeKutta4Step(BODY *body, CONTROL *control, SYSTEM *system,
2 * evolve->daDeriv[2][iBody][iVar] +
evolve->daDeriv[3][iBody][iVar]);
for (iEqn = 0; iEqn < update[iBody].iNumEqns[iVar]; iEqn++) {
update[iBody].daDerivProc[iVar][iEqn] = 1. / 6 *
(evolve->daDerivProc[0][iBody][iVar][iEqn] +
2 * evolve->daDerivProc[1][iBody][iVar][iEqn] +
2 * evolve->daDerivProc[2][iBody][iVar][iEqn] +
evolve->daDerivProc[3][iBody][iVar][iEqn]);
update[iBody].daDerivProc[iVar][iEqn] =
1. / 6 *
(evolve->daDerivProc[0][iBody][iVar][iEqn] +
2 * evolve->daDerivProc[1][iBody][iVar][iEqn] +
2 * evolve->daDerivProc[2][iBody][iVar][iEqn] +
evolve->daDerivProc[3][iBody][iVar][iEqn]);
}

if (update[iBody].iaType[iVar][0] == 0 ||
Expand All @@ -604,8 +609,8 @@ void Evolve(BODY *body, CONTROL *control, FILES *files, MODULE *module,
fnIntegrate fnOneStep) {
/* Master evolution routine that controls the simulation integration. */
int iDir, iBody, iModule; // Dummy counting variables
double dDt, dFoo; // Next timestep, dummy variable
double dEqSpinRate; // Store the equilibrium spin rate
double dDt, dFoo; // Next timestep, dummy variable
double dEqSpinRate; // Store the equilibrium spin rate

if (control->Evolve.bDoForward) {
iDir = 1;
Expand Down Expand Up @@ -644,7 +649,7 @@ void Evolve(BODY *body, CONTROL *control, FILES *files, MODULE *module,
*/

control->Evolve.iStepsSinceLastOutput = 0;
control->Evolve.iTotalSteps = 0;
control->Evolve.iTotalSteps = 0;
while (control->Evolve.dTime < control->Evolve.dStopTime) {
/* Take one step */
fnOneStep(body, control, system, update, fnUpdate, &dDt, iDir);
Expand Down Expand Up @@ -687,7 +692,10 @@ void Evolve(BODY *body, CONTROL *control, FILES *files, MODULE *module,
// Timesteps are synchronized with the output time, so this statement is
// sufficient
control->Io.dNextOutput += control->Io.dOutputTime;
control->Evolve.iStepsSinceLastOutput = 0;
if (control->Evolve.dTime < control->Evolve.dStopTime) {
// Reset counter if we are not at the end of the simulation
control->Evolve.iStepsSinceLastOutput = 0;
}
}

/* Get auxiliary properties for next step -- first call
Expand Down
Loading

0 comments on commit d15e244

Please sign in to comment.