Skip to content
Open
Show file tree
Hide file tree
Changes from 5 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
28 changes: 14 additions & 14 deletions libLanlGeoMag/ComputeI_FromMltMlat.c
Original file line number Diff line number Diff line change
Expand Up @@ -34,7 +34,7 @@ double ComputeI_FromMltMlat1( double Bm, double MLT, double mlat, double *r, dou

int reset=1, reset2;

double I, Phi, cl, sl, rat, SS1, SS2, SS, Sn, Ss, Htry, Hdid, Hnext, Bs, Be, s, sgn;
double I1, Phi, cl, sl, rat, SS1, SS2, SS, Sn, Ss, Htry, Hdid, Hnext, Bs, Be, s, sgn;
Lgm_Vector w, u, Pmirror1, Pmirror2, v1, v2, v3, Bvec, P, Ps, u_scale, Bvectmp, Ptmp;
double stmp, Btmp;

Expand All @@ -52,8 +52,8 @@ double ComputeI_FromMltMlat1( double Bm, double MLT, double mlat, double *r, dou
* Couldnt get a valid Bm. (The bracket is pretty huge,so
* we probably ought to believe there really isnt a valid one.)
*/
if (LstarInfo->VerbosityLevel > 1) printf("\t%sNo Bm found: setting I to 9e99%s\n", LstarInfo->PreStr, LstarInfo->PostStr);
I = 9e99;
if (LstarInfo->VerbosityLevel > 1) printf("\t%sNo Bm found: setting I1 to 9e99%s\n", LstarInfo->PreStr, LstarInfo->PostStr);
I1 = 9e99;

} else {

Expand Down Expand Up @@ -195,7 +195,7 @@ double ComputeI_FromMltMlat1( double Bm, double MLT, double mlat, double *r, dou
*
* Trace from Pm_North to Pm_South
*/
I = 9e99;
I1 = 9e99;
//LstarInfo->mInfo->Hmax = 10.0;
//LstarInfo->mInfo->Hmax = 0.1;
LstarInfo->mInfo->Hmax = 0.1;
Expand Down Expand Up @@ -306,18 +306,18 @@ if (0==1){
/*
* Do I integral with interped integrand.
*/
//printf("I = %g\n", I);
I = Iinv_interped( LstarInfo->mInfo );
//printf("I = %g Sm_South, Sm_North = %g %g\n", I, LstarInfo->mInfo->Sm_South, LstarInfo->mInfo->Sm_North);
// if (LstarInfo->VerbosityLevel > 1) printf("\t\t%s Integral Invariant, I (interped): %15.8g I-I0: %15.8g [a,b]: %.15g %.15g mlat: %12.8lf (nCalls = %d)%s\n", LstarInfo->PreStr, I, I-I0, LstarInfo->mInfo->Sm_South, LstarInfo->mInfo->Sm_North, mlat, LstarInfo->mInfo->Lgm_n_I_integrand_Calls, LstarInfo->PostStr );
//printf("I = %g\n", I1);
I1 = Iinv_interped( LstarInfo->mInfo );
Comment thread
drsteve marked this conversation as resolved.
//printf("I = %g Sm_South, Sm_North = %g %g\n", I1, LstarInfo->mInfo->Sm_South, LstarInfo->mInfo->Sm_North);
// if (LstarInfo->VerbosityLevel > 1) printf("\t\t%s Integral Invariant, I (interped): %15.8g I-I0: %15.8g [a,b]: %.15g %.15g mlat: %12.8lf (nCalls = %d)%s\n", LstarInfo->PreStr, I1, I1-I0, LstarInfo->mInfo->Sm_South, LstarInfo->mInfo->Sm_North, mlat, LstarInfo->mInfo->Lgm_n_I_integrand_Calls, LstarInfo->PostStr );
if (LstarInfo->VerbosityLevel > 1) {
printf("\t\t%s mlat: %13.6g I: %13.6g I0: %13.6g I-I0: %13.6g [Sa,Sb]: %.8g %.8g (nCalls = %d)%s\n", LstarInfo->PreStr, mlat, I, I0, I-I0, LstarInfo->mInfo->Sm_South, LstarInfo->mInfo->Sm_North, LstarInfo->mInfo->Lgm_n_I_integrand_Calls, LstarInfo->PostStr );
printf("\t\t%s mlat: %13.6g I1: %13.6g I0: %13.6g I1-I0: %13.6g [Sa,Sb]: %.8g %.8g (nCalls = %d)%s\n", LstarInfo->PreStr, mlat, I1, I0, I1-I0, LstarInfo->mInfo->Sm_South, LstarInfo->mInfo->Sm_North, LstarInfo->mInfo->Lgm_n_I_integrand_Calls, LstarInfo->PostStr );
}
FreeSpline( LstarInfo->mInfo );

} else {

I = 9e99;
I1 = 9e99;

}

Expand All @@ -332,19 +332,19 @@ if (0==1){
/*
* Do full blown I integral.
*/
I = Iinv( LstarInfo->mInfo );
if (LstarInfo->VerbosityLevel > 1) printf("\t\t%s Integral Invariant, I (full integral): %15.8g I-I0: %15.8g mlat: %12.8lf (nCalls = %d)%s\n", LstarInfo->PreStr, I, I-I0, mlat, LstarInfo->mInfo->Lgm_n_I_integrand_Calls, LstarInfo->PostStr );
I1 = Iinv( LstarInfo->mInfo );
if (LstarInfo->VerbosityLevel > 1) printf("\t\t%s Integral Invariant, I1 (full integral): %15.8g I1-I0: %15.8g mlat: %12.8lf (nCalls = %d)%s\n", LstarInfo->PreStr, I1, I1-I0, mlat, LstarInfo->mInfo->Lgm_n_I_integrand_Calls, LstarInfo->PostStr );
}

} else {
I = 9e99;
I1 = 9e99;
}

}




return( I );
return( I1 );

}
40 changes: 20 additions & 20 deletions libLanlGeoMag/ComputeI_FromMltMlat2.c
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@ double ComputeI_FromMltMlat2( double Bm, double MLT, double mlat, double *r, dou

int reset=1, reset2, TraceFlag;

double Bmin, I, Phi, cl, sl, rat, SS1, SS2, SS, Sn, Ss, Htry, Hdid, Hnext, Bs, Be, s, sgn;
double Bmin, I1, Phi, cl, sl, rat, SS1, SS2, SS, Sn, Ss, Htry, Hdid, Hnext, Bs, Be, s, sgn;
Lgm_Vector w, u, Pmirror1, Pmirror2, v1, v2, v3, Bvec, P, Ps, u_scale, Bvectmp, Ptmp;
double stmp, Btmp;

Expand Down Expand Up @@ -38,13 +38,13 @@ double ComputeI_FromMltMlat2( double Bm, double MLT, double mlat, double *r, dou

if ( TraceFlag != LGM_CLOSED ) {

I = 9e99;
I1 = 9e99;
if (LstarInfo->VerbosityLevel > 1) {
printf("\t\t%s mlat: %13.6g I: %13.6g I0: %13.6g > Field Line not closed%s\n", LstarInfo->PreStr, mlat, I, I0, LstarInfo->PostStr );
printf("\t\t%s mlat: %13.6g I1: %13.6g I0: %13.6g > Field Line not closed%s\n", LstarInfo->PreStr, mlat, I1, I0, LstarInfo->PostStr );
}

*ErrorStatus = -1; // FL open. I undefined.
return( I );
return( I1 );

} else if ( Bmin <= Bm ) {

Expand All @@ -67,40 +67,40 @@ double ComputeI_FromMltMlat2( double Bm, double MLT, double mlat, double *r, dou
//printf("Pmirror2 = %g %g %g\n", Pmirror2.x, Pmirror2.y, Pmirror2.z);

} else {
I = 9e99;
I1 = 9e99;
if (LstarInfo->VerbosityLevel > 1) {
printf("\t\t%s mlat: %13.6g I: %13.6g I0: %13.6g > Unable to find southern mirror point%s\n", LstarInfo->PreStr, mlat, I, I0, LstarInfo->PostStr );
printf("\t\t%s mlat: %13.6g I1: %13.6g I0: %13.6g > Unable to find southern mirror point%s\n", LstarInfo->PreStr, mlat, I1, I0, LstarInfo->PostStr );
}
*ErrorStatus = -3; // No valid mirror point in the south
return( I );
return( I1 );
}

// total distance between mirror point.
SS = SS1 + SS2;

// If its really small, just return 0.0 for I
if ( fabs(SS) < 1e-7) {
I = 0.0;
I1 = 0.0;
if (LstarInfo->VerbosityLevel > 1) {
printf("\t\t%s mlat: %13.6g I: %13.6g I0: %13.6g > Distance between mirror points is < 1e-7Re, Assuming I=0%s\n", LstarInfo->PreStr, mlat, I, I0, LstarInfo->PostStr );
printf("\t\t%s mlat: %13.6g I1: %13.6g I0: %13.6g > Distance between mirror points is < 1e-7Re, Assuming I1=0%s\n", LstarInfo->PreStr, mlat, I1, I0, LstarInfo->PostStr );
}
*ErrorStatus = 2; // Flag that we did this
return( I );
return( I1 );
}

} else {
I = 9e99;
I1 = 9e99;
if (LstarInfo->VerbosityLevel > 1) {
printf("\t\t%s mlat: %13.6g I: %13.6g I0: %13.6g > Unable to find northern mirror point%s\n", LstarInfo->PreStr, mlat, I, I0, LstarInfo->PostStr );
printf("\t\t%s mlat: %13.6g I1: %13.6g I0: %13.6g > Unable to find northern mirror point%s\n", LstarInfo->PreStr, mlat, I1, I0, LstarInfo->PostStr );
}
*ErrorStatus = -2; // No valid mirror point in the north
return( I );
return( I1 );
}

/*
* OK, we have both mirror points. Lets compute I
*/
I = 9e99;
I1 = 9e99;
LstarInfo->mInfo->Hmax = 0.1;
LstarInfo->mInfo->Hmax = SS/(double)LstarInfo->mInfo->nDivs;

Expand Down Expand Up @@ -137,16 +137,16 @@ double ComputeI_FromMltMlat2( double Bm, double MLT, double mlat, double *r, dou
/*
* Do I integral with interped integrand.
*/
I = Iinv_interped( LstarInfo->mInfo );
I1 = Iinv_interped( LstarInfo->mInfo );
if (LstarInfo->VerbosityLevel > 1) {
printf("\t\t%s mlat: %13.6g I: %13.6g I0: %13.6g I-I0: %13.6g [Sa,Sb]: %.8g %.8g (nCalls = %d)%s\n", LstarInfo->PreStr, mlat, I, I0, I-I0, LstarInfo->mInfo->Sm_South, LstarInfo->mInfo->Sm_North, LstarInfo->mInfo->Lgm_n_I_integrand_Calls, LstarInfo->PostStr );
printf("\t\t%s mlat: %13.6g I1: %13.6g I0: %13.6g I1-I0: %13.6g [Sa,Sb]: %.8g %.8g (nCalls = %d)%s\n", LstarInfo->PreStr, mlat, I1, I0, I1-I0, LstarInfo->mInfo->Sm_South, LstarInfo->mInfo->Sm_North, LstarInfo->mInfo->Lgm_n_I_integrand_Calls, LstarInfo->PostStr );
}
FreeSpline( LstarInfo->mInfo );

} else {
I = 9e99;
I1 = 9e99;
if (LstarInfo->VerbosityLevel > 1) {
printf("\t\t%s mlat: %13.6g I: %13.6g I0: %13.6g Couldnt initialize spline%s\n", LstarInfo->PreStr, mlat, I, I0, LstarInfo->PostStr );
printf("\t\t%s mlat: %13.6g I1: %13.6g I0: %13.6g Couldnt initialize spline%s\n", LstarInfo->PreStr, mlat, I1, I0, LstarInfo->PostStr );
}
}

Expand All @@ -160,14 +160,14 @@ double ComputeI_FromMltMlat2( double Bm, double MLT, double mlat, double *r, dou

*ErrorStatus = -10; // We didnt even try to compute I because we Field line min-B is greater than Bm.
if (LstarInfo->VerbosityLevel > 1) {
printf("\t\t%s mlat: %13.6g I: undefined I0: %13.6g I undefined. min-B is greater than Bm.%s\n", LstarInfo->PreStr, mlat, I0, LstarInfo->PostStr );
printf("\t\t%s mlat: %13.6g I1: undefined I0: %13.6g I1 undefined. min-B is greater than Bm.%s\n", LstarInfo->PreStr, mlat, I0, LstarInfo->PostStr );
}
return( 9e99 );

}


return( I );
return( I1 );


}
Loading