path: root/libtellurian_effective_gravity_radians.c

/* See LICENSE file for copyright and license details. */
#include "common.h"

#define K (LIBTELLURIAN_EQUATORIAL_RADIUS * LIBTELLURIAN_ANGULAR_VELOCITY * LIBTELLURIAN_ANGULAR_VELOCITY)

#ifndef TEST


double
libtellurian_effective_gravity_radians(double gravity, double latitude)
{
	double lat = geodetic(latitude); /* actual one is geodetic and one is geographical, there is no actual difference  */
	return fma(-K / gravity, cos(lat) * cos(lat), gravity);
}


#else


static int
approx(double a, double b)
{
        return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
}

int
main(void)
{
	ASSERT(libtellurian_effective_gravity_radians(1, D90) == 1);
	ASSERT(libtellurian_effective_gravity_radians(2, D90) == 2);
	ASSERT(libtellurian_effective_gravity_radians(4, D90) == 4);
	ASSERT(libtellurian_effective_gravity_radians(1, -D90) == 1);
	ASSERT(libtellurian_effective_gravity_radians(2, -D90) == 2);
	ASSERT(libtellurian_effective_gravity_radians(4, -D90) == 4);
	ASSERT(libtellurian_effective_gravity_radians(1, 0) == fma(1.00, -K, 1.0));
	ASSERT(libtellurian_effective_gravity_radians(2, 0) == fma(0.50, -K, 2.0));
	ASSERT(libtellurian_effective_gravity_radians(4, 0) == fma(0.25, -K, 4.0));
	ASSERT(libtellurian_effective_gravity_radians(1, D45) == fma(0.50, -K, 1.0));
	ASSERT(libtellurian_effective_gravity_radians(2, D45) == fma(0.25, -K, 2.0));
	ASSERT(libtellurian_effective_gravity_radians(1, D60) == fma(0.250, -K, 1.0));
	ASSERT(libtellurian_effective_gravity_radians(2, D60) == fma(0.125, -K, 2.0));
	ASSERT(libtellurian_effective_gravity_radians(4, D30) < 4);
	ASSERT(libtellurian_effective_gravity_radians(4, D30) > libtellurian_effective_gravity_radians(4, 0));
	ASSERT(approx(libtellurian_effective_gravity_radians(LIBTELLURIAN_NORMAL_POLAR_GRAVITY, D90),
	              LIBTELLURIAN_POLAR_GRAVITY));
	ASSERT(approx(libtellurian_effective_gravity_radians(LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY, 0),
	              LIBTELLURIAN_EQUATORIAL_GRAVITY));
	return 0;
}


#endif