Jupiter Moons: stubbed the idea of a model optimizer.

I want to experiment with truncating the L1.2 series to
sacrifice some accuracy for smaller generated code.
To that end, I implemented the ability to save the
Jupiter moons model after loading it. I added a 'jmopt'
command to the 'generate' program that will do this
optimization. For now, it just loads the model and
saves it back to a different file. Then the code generator
loads from the saved file instead of the original.
This commit verifies that everything is still working,
before I start truncating the series.
This commit is contained in:
Don Cross
2021-04-11 13:42:23 -04:00
parent 87f5dd1c95
commit 43d05b25bd
5 changed files with 140 additions and 35 deletions

View File

@@ -28,7 +28,6 @@
#include <stdarg.h>
#include "novas.h"
#include "codegen.h"
#include "vsop.h"
#include "top2013.h"
#include "ephfile.h"
@@ -1380,38 +1379,6 @@ fail:
/*-------------------------- begin Jupiter moons -----------------------------*/
#define NUM_JUPITER_MOONS 4
#define MAX_JM_SERIES 50
#define NUM_JM_VARS 4
#define MAX_JUPITER_TERMS (NUM_JUPITER_MOONS * NUM_JM_VARS * MAX_JM_SERIES)
typedef struct
{
double mu; /* mu = G(M+m), where M = Jupiter mass, m = moon mass. */
double al[2]; /* mean longitude coefficients */
vsop_series_t a;
vsop_series_t l;
vsop_series_t z;
vsop_series_t zeta;
}
jupiter_moon_t;
typedef struct
{
/* angles that rotate Jupiter equatorial system to EQJ */
double incl;
double psi;
/* a rotation matrix calculated from incl and psi, for the convenience of the code generators. */
double rot[3][3];
jupiter_moon_t moon[NUM_JUPITER_MOONS];
vsop_term_t buffer[MAX_JUPITER_TERMS];
}
jupiter_moon_model_t;
static int ReadFixLine(char *line, size_t size, FILE *infile)
{
if (!fgets(line, size, infile))
@@ -1433,7 +1400,7 @@ static int ReadFixLine(char *line, size_t size, FILE *infile)
do{if ((nscanned) != (required)) FAIL("LoadJupiterMoonModel(%s line %d): expected %d tokens, found %d\n", filename, lnum, required, nscanned);}while(0)
static int LoadJupiterMoonModel(const char *filename, jupiter_moon_model_t *model)
int LoadJupiterMoonModel(const char *filename, jupiter_moon_model_t *model)
{
int error;
FILE *infile;
@@ -1572,6 +1539,65 @@ fail:
}
static int JupiterMoonWriteVar(FILE *outfile, const vsop_series_t *series, const double *al)
{
int tindex;
fprintf(outfile, "%d\n", series->nterms_calc);
if (al != NULL)
fprintf(outfile, "%23.16le %23.16le\n", al[0], al[1]);
for (tindex = 0; tindex < series->nterms_calc; ++tindex)
{
fprintf(outfile, "%3d %23.16lf %20.13le %20.13le\n",
1 + tindex,
series->term[tindex].amplitude,
series->term[tindex].phase,
series->term[tindex].frequency);
}
return 0;
}
int SaveJupiterMoonModel(const char *filename, const jupiter_moon_model_t *model)
{
int error, mindex;
FILE *outfile;
outfile = fopen(filename, "wt");
if (outfile == NULL)
FAIL("SaveJupiterMoonModel: cannot open output file: %s\n", filename);
/* As a hack, we emit 18 blank lines, because we ignore those from the original L1.2 file. */
for (mindex = 0; mindex < 18; ++mindex)
fprintf(outfile, "\n");
/* Write G(M+m) constants */
fprintf(outfile, "%23.16le %23.16le %23.16le %23.16le\n", model->moon[0].mu, model->moon[1].mu, model->moon[2].mu, model->moon[3].mu);
/* Write Jupiter equatorial orientation angles. */
fprintf(outfile, "%23.16le %23.16le\n", model->psi, model->incl);
for (mindex = 0; mindex < NUM_JUPITER_MOONS; ++mindex)
{
CHECK(JupiterMoonWriteVar(outfile, &model->moon[mindex].a, NULL));
CHECK(JupiterMoonWriteVar(outfile, &model->moon[mindex].l, model->moon[mindex].al));
CHECK(JupiterMoonWriteVar(outfile, &model->moon[mindex].z, NULL));
CHECK(JupiterMoonWriteVar(outfile, &model->moon[mindex].zeta, NULL));
}
fprintf(outfile, "\n"); /* final blank line needed for loader to know it has all 4 moons. */
error = 0;
fail:
if (outfile != NULL) fclose(outfile);
return error;
}
static int JupiterMoons_C(cg_context_t *context, const jupiter_moon_model_t *model)
{
int mindex, var, i;
@@ -1658,7 +1684,7 @@ static int JupiterMoons(cg_context_t *context)
if (model == NULL)
FAIL("JupiterMoons: memory allocation failure!\n");
CHECK(LoadJupiterMoonModel("jupiter_moons/fortran/BisL1.2.dat", model));
CHECK(LoadJupiterMoonModel("output/jupiter_moons.txt", model));
switch (context->language)
{

View File

@@ -24,6 +24,8 @@
#ifndef __DDC_ASTRO_CODEGEN_H
#define __DDC_ASTRO_CODEGEN_H
#include "vsop.h"
#define CHECK(x) do{if(0 != (error = (x))) goto fail;}while(0)
#define FAIL(...) do{fprintf(stderr, __VA_ARGS__); error = 1; goto fail;}while(0)
@@ -45,4 +47,41 @@ int GenerateCode(
double ExtrapolatedDeltaT(int year);
/*-------------------------- Jupiter moons -----------------------------*/
#define NUM_JUPITER_MOONS 4
#define MAX_JM_SERIES 50
#define NUM_JM_VARS 4
#define MAX_JUPITER_TERMS (NUM_JUPITER_MOONS * NUM_JM_VARS * MAX_JM_SERIES)
typedef struct
{
double mu; /* mu = G(M+m), where M = Jupiter mass, m = moon mass. */
double al[2]; /* mean longitude coefficients */
vsop_series_t a;
vsop_series_t l;
vsop_series_t z;
vsop_series_t zeta;
}
jupiter_moon_t;
typedef struct
{
/* angles that rotate Jupiter equatorial system to EQJ */
double incl;
double psi;
/* a rotation matrix calculated from incl and psi, for the convenience of the code generators. */
double rot[3][3];
jupiter_moon_t moon[NUM_JUPITER_MOONS];
vsop_term_t buffer[MAX_JUPITER_TERMS];
}
jupiter_moon_model_t;
int LoadJupiterMoonModel(const char *filename, jupiter_moon_model_t *model);
int SaveJupiterMoonModel(const char *filename, const jupiter_moon_model_t *model);
#endif /* __DDC_ASTRO_CODEGEN_H */

View File

@@ -79,6 +79,7 @@ static int PrintUsage(void);
static int GenerateVsopPlanets(void);
static int TopCalc(const char *name, const char *date);
static int GenerateApsisTestData(void);
static int OptimizeJupiterMoons(const char *inFileName, const char *outFileName);
static int GenerateSource(void);
static int TestVsopModel(vsop_model_t *model, int body, double threshold, double *max_arcmin, int *trunc_terms);
static int SaveVsopFile(const vsop_model_t *model);
@@ -144,6 +145,9 @@ int main(int argc, const char *argv[])
if (argc == 2 && !strcmp(argv[1], "apsis"))
return GenerateApsisTestData();
if (argc == 2 && !strcmp(argv[1], "jmopt"))
return OptimizeJupiterMoons("jupiter_moons/fortran/BisL1.2.dat", "output/jupiter_moons.txt");
if (argc == 2 && !strcmp(argv[1], "source"))
return GenerateSource();
@@ -1849,3 +1853,25 @@ fail:
return error;
}
/*------------------------------------------------------------------------------------------------*/
static int OptimizeJupiterMoons(const char *inFileName, const char *outFileName)
{
int error;
jupiter_moon_model_t *model;
model = calloc(1, sizeof(jupiter_moon_model_t));
if (model == NULL)
FAIL("OptimizeJupiterMoons: memory allocation failure.\n");
CHECK(LoadJupiterMoonModel(inFileName, model));
CHECK(SaveJupiterMoonModel(outFileName, model));
error = 0;
fail:
free(model);
return error;
}
/*------------------------------------------------------------------------------------------------*/

View File

@@ -107,6 +107,11 @@ if [[ "${FASTMODE}" == "false" ]]; then
./generate planets || Fail "Could not generate planet models."
fi
echo ""
echo "Optimizing Jupiter Moon models."
rm -f output/jupiter_moons.txt
./generate jmopt || Fail "Error optimizing Jupiter Moon models."
echo ""
echo "Generating apsis test data."
rm -f apsides/apsis_*.txt

View File

@@ -78,6 +78,15 @@ if !FASTMODE! == false (
)
)
echo.
echo.Optimizing Jupiter Moon models.
if exist output\jupiter_moons.txt (del output\jupiter_moons.txt)
!GENEXE! jmopt
if errorlevel 1 (
echo.FATAL: !GENEXE! jmopt
exit /b 1
)
echo.
echo.Generating apsis test data.
if exist apsides\apsis_*.txt (del apsides\apsis_*.txt)