#include <tgmath.h>
#include <gmp.h>
#include <mpfr.h>
#include <stdio.h>
#include <stdbool.h>

char *name = "Henrywood and Agarwal, Equation (9a)";

double f_if(float w0, float M, float D, float h, float l, float d) {
        float r31953 = w0;
        float r31954 = 1;
        float r31955 = M;
        float r31956 = D;
        float r31957 = r31955 * r31956;
        float r31958 = 2;
        float r31959 = d;
        float r31960 = r31958 * r31959;
        float r31961 = r31957 / r31960;
        float r31962 = pow(r31961, r31958);
        float r31963 = h;
        float r31964 = l;
        float r31965 = r31963 / r31964;
        float r31966 = r31962 * r31965;
        float r31967 = r31954 - r31966;
        float r31968 = sqrt(r31967);
        float r31969 = r31953 * r31968;
        return r31969;
}

double f_id(double w0, double M, double D, double h, double l, double d) {
        double r31970 = w0;
        double r31971 = 1;
        double r31972 = M;
        double r31973 = D;
        double r31974 = r31972 * r31973;
        double r31975 = 2;
        double r31976 = d;
        double r31977 = r31975 * r31976;
        double r31978 = r31974 / r31977;
        double r31979 = pow(r31978, r31975);
        double r31980 = h;
        double r31981 = l;
        double r31982 = r31980 / r31981;
        double r31983 = r31979 * r31982;
        double r31984 = r31971 - r31983;
        double r31985 = sqrt(r31984);
        double r31986 = r31970 * r31985;
        return r31986;
}


double f_of(float w0, float M, float D, float h, float l, float d) {
        float r31987 = w0;
        float r31988 = 1;
        float r31989 = M;
        float r31990 = D;
        float r31991 = r31989 * r31990;
        float r31992 = 2;
        float r31993 = d;
        float r31994 = r31992 * r31993;
        float r31995 = r31991 / r31994;
        float r31996 = h;
        float r31997 = r31995 * r31996;
        float r31998 = r31995 * r31997;
        float r31999 = l;
        float r32000 = r31998 / r31999;
        float r32001 = r31988 - r32000;
        float r32002 = sqrt(r32001);
        float r32003 = r31987 * r32002;
        return r32003;
}

double f_od(double w0, double M, double D, double h, double l, double d) {
        double r32004 = w0;
        double r32005 = 1;
        double r32006 = M;
        double r32007 = D;
        double r32008 = r32006 * r32007;
        double r32009 = 2;
        double r32010 = d;
        double r32011 = r32009 * r32010;
        double r32012 = r32008 / r32011;
        double r32013 = h;
        double r32014 = r32012 * r32013;
        double r32015 = r32012 * r32014;
        double r32016 = l;
        double r32017 = r32015 / r32016;
        double r32018 = r32005 - r32017;
        double r32019 = sqrt(r32018);
        double r32020 = r32004 * r32019;
        return r32020;
}

void mpfr_fmod2(mpfr_t r, mpfr_t n, mpfr_t d, mpfr_rnd_t rmd) {
        mpfr_fmod(r, n, d, rmd);
        if (mpfr_cmp_ui(r, 0) < 0) mpfr_add(r, r, d, rmd);
}


static mpfr_t r32021, r32022, r32023, r32024, r32025, r32026, r32027, r32028, r32029, r32030, r32031, r32032, r32033, r32034, r32035, r32036, r32037;

void setup_mpfr_f_im() {
        mpfr_set_default_prec(400);
        mpfr_init(r32021);
        mpfr_init_set_str(r32022, "1", 10, MPFR_RNDN);
        mpfr_init(r32023);
        mpfr_init(r32024);
        mpfr_init(r32025);
        mpfr_init_set_str(r32026, "2", 10, MPFR_RNDN);
        mpfr_init(r32027);
        mpfr_init(r32028);
        mpfr_init(r32029);
        mpfr_init(r32030);
        mpfr_init(r32031);
        mpfr_init(r32032);
        mpfr_init(r32033);
        mpfr_init(r32034);
        mpfr_init(r32035);
        mpfr_init(r32036);
        mpfr_init(r32037);
}

double f_im(double w0, double M, double D, double h, double l, double d) {
        mpfr_set_d(r32021, w0, MPFR_RNDN);
        ;
        mpfr_set_d(r32023, M, MPFR_RNDN);
        mpfr_set_d(r32024, D, MPFR_RNDN);
        mpfr_mul(r32025, r32023, r32024, MPFR_RNDN);
        ;
        mpfr_set_d(r32027, d, MPFR_RNDN);
        mpfr_mul(r32028, r32026, r32027, MPFR_RNDN);
        mpfr_div(r32029, r32025, r32028, MPFR_RNDN);
        mpfr_pow(r32030, r32029, r32026, MPFR_RNDN);
        mpfr_set_d(r32031, h, MPFR_RNDN);
        mpfr_set_d(r32032, l, MPFR_RNDN);
        mpfr_div(r32033, r32031, r32032, MPFR_RNDN);
        mpfr_mul(r32034, r32030, r32033, MPFR_RNDN);
        mpfr_sub(r32035, r32022, r32034, MPFR_RNDN);
        mpfr_sqrt(r32036, r32035, MPFR_RNDN);
        mpfr_mul(r32037, r32021, r32036, MPFR_RNDN);
        return mpfr_get_d(r32037, MPFR_RNDN);
}

static mpfr_t r32038, r32039, r32040, r32041, r32042, r32043, r32044, r32045, r32046, r32047, r32048, r32049, r32050, r32051, r32052, r32053, r32054;

void setup_mpfr_f_fm() {
        mpfr_set_default_prec(400);
        mpfr_init(r32038);
        mpfr_init_set_str(r32039, "1", 10, MPFR_RNDN);
        mpfr_init(r32040);
        mpfr_init(r32041);
        mpfr_init(r32042);
        mpfr_init_set_str(r32043, "2", 10, MPFR_RNDN);
        mpfr_init(r32044);
        mpfr_init(r32045);
        mpfr_init(r32046);
        mpfr_init(r32047);
        mpfr_init(r32048);
        mpfr_init(r32049);
        mpfr_init(r32050);
        mpfr_init(r32051);
        mpfr_init(r32052);
        mpfr_init(r32053);
        mpfr_init(r32054);
}

double f_fm(double w0, double M, double D, double h, double l, double d) {
        mpfr_set_d(r32038, w0, MPFR_RNDN);
        ;
        mpfr_set_d(r32040, M, MPFR_RNDN);
        mpfr_set_d(r32041, D, MPFR_RNDN);
        mpfr_mul(r32042, r32040, r32041, MPFR_RNDN);
        ;
        mpfr_set_d(r32044, d, MPFR_RNDN);
        mpfr_mul(r32045, r32043, r32044, MPFR_RNDN);
        mpfr_div(r32046, r32042, r32045, MPFR_RNDN);
        mpfr_set_d(r32047, h, MPFR_RNDN);
        mpfr_mul(r32048, r32046, r32047, MPFR_RNDN);
        mpfr_mul(r32049, r32046, r32048, MPFR_RNDN);
        mpfr_set_d(r32050, l, MPFR_RNDN);
        mpfr_div(r32051, r32049, r32050, MPFR_RNDN);
        mpfr_sub(r32052, r32039, r32051, MPFR_RNDN);
        mpfr_sqrt(r32053, r32052, MPFR_RNDN);
        mpfr_mul(r32054, r32038, r32053, MPFR_RNDN);
        return mpfr_get_d(r32054, MPFR_RNDN);
}

static mpfr_t r32055, r32056, r32057, r32058, r32059, r32060, r32061, r32062, r32063, r32064, r32065, r32066, r32067, r32068, r32069, r32070, r32071;

void setup_mpfr_f_dm() {
        mpfr_set_default_prec(400);
        mpfr_init(r32055);
        mpfr_init_set_str(r32056, "1", 10, MPFR_RNDN);
        mpfr_init(r32057);
        mpfr_init(r32058);
        mpfr_init(r32059);
        mpfr_init_set_str(r32060, "2", 10, MPFR_RNDN);
        mpfr_init(r32061);
        mpfr_init(r32062);
        mpfr_init(r32063);
        mpfr_init(r32064);
        mpfr_init(r32065);
        mpfr_init(r32066);
        mpfr_init(r32067);
        mpfr_init(r32068);
        mpfr_init(r32069);
        mpfr_init(r32070);
        mpfr_init(r32071);
}

double f_dm(double w0, double M, double D, double h, double l, double d) {
        mpfr_set_d(r32055, w0, MPFR_RNDN);
        ;
        mpfr_set_d(r32057, M, MPFR_RNDN);
        mpfr_set_d(r32058, D, MPFR_RNDN);
        mpfr_mul(r32059, r32057, r32058, MPFR_RNDN);
        ;
        mpfr_set_d(r32061, d, MPFR_RNDN);
        mpfr_mul(r32062, r32060, r32061, MPFR_RNDN);
        mpfr_div(r32063, r32059, r32062, MPFR_RNDN);
        mpfr_set_d(r32064, h, MPFR_RNDN);
        mpfr_mul(r32065, r32063, r32064, MPFR_RNDN);
        mpfr_mul(r32066, r32063, r32065, MPFR_RNDN);
        mpfr_set_d(r32067, l, MPFR_RNDN);
        mpfr_div(r32068, r32066, r32067, MPFR_RNDN);
        mpfr_sub(r32069, r32056, r32068, MPFR_RNDN);
        mpfr_sqrt(r32070, r32069, MPFR_RNDN);
        mpfr_mul(r32071, r32055, r32070, MPFR_RNDN);
        return mpfr_get_d(r32071, MPFR_RNDN);
}

