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

char *name = "Jmat.Real.lambertw, newton loop step";

double f_if(float wj, float x) {
        float r28040 = wj;
        float r28041 = exp(r28040);
        float r28042 = r28040 * r28041;
        float r28043 = x;
        float r28044 = r28042 - r28043;
        float r28045 = r28041 + r28042;
        float r28046 = r28044 / r28045;
        float r28047 = r28040 - r28046;
        return r28047;
}

double f_id(double wj, double x) {
        double r28048 = wj;
        double r28049 = exp(r28048);
        double r28050 = r28048 * r28049;
        double r28051 = x;
        double r28052 = r28050 - r28051;
        double r28053 = r28049 + r28050;
        double r28054 = r28052 / r28053;
        double r28055 = r28048 - r28054;
        return r28055;
}


double f_of(float wj, float x) {
        float r28056 = wj;
        float r28057 = r28056 * r28056;
        float r28058 = r28056 - r28057;
        float r28059 = 1;
        float r28060 = r28059 / r28056;
        float r28061 = -4;
        float r28062 = pow(r28060, r28061);
        float r28063 = fma(r28056, r28058, r28062);
        float r28064 = x;
        float r28065 = exp(r28056);
        float r28066 = fma(r28056, r28065, r28065);
        float r28067 = r28064 / r28066;
        float r28068 = r28063 + r28067;
        return r28068;
}

double f_od(double wj, double x) {
        double r28069 = wj;
        double r28070 = r28069 * r28069;
        double r28071 = r28069 - r28070;
        double r28072 = 1;
        double r28073 = r28072 / r28069;
        double r28074 = -4;
        double r28075 = pow(r28073, r28074);
        double r28076 = fma(r28069, r28071, r28075);
        double r28077 = x;
        double r28078 = exp(r28069);
        double r28079 = fma(r28069, r28078, r28078);
        double r28080 = r28077 / r28079;
        double r28081 = r28076 + r28080;
        return r28081;
}

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 r28082, r28083, r28084, r28085, r28086, r28087, r28088, r28089;

void setup_mpfr_f_im() {
        mpfr_set_default_prec(848);
        mpfr_init(r28082);
        mpfr_init(r28083);
        mpfr_init(r28084);
        mpfr_init(r28085);
        mpfr_init(r28086);
        mpfr_init(r28087);
        mpfr_init(r28088);
        mpfr_init(r28089);
}

double f_im(double wj, double x) {
        mpfr_set_d(r28082, wj, MPFR_RNDN);
        mpfr_exp(r28083, r28082, MPFR_RNDN);
        mpfr_mul(r28084, r28082, r28083, MPFR_RNDN);
        mpfr_set_d(r28085, x, MPFR_RNDN);
        mpfr_sub(r28086, r28084, r28085, MPFR_RNDN);
        mpfr_add(r28087, r28083, r28084, MPFR_RNDN);
        mpfr_div(r28088, r28086, r28087, MPFR_RNDN);
        mpfr_sub(r28089, r28082, r28088, MPFR_RNDN);
        return mpfr_get_d(r28089, MPFR_RNDN);
}

static mpfr_t r28090, r28091, r28092, r28093, r28094, r28095, r28096, r28097, r28098, r28099, r28100, r28101, r28102;

void setup_mpfr_f_fm() {
        mpfr_set_default_prec(848);
        mpfr_init(r28090);
        mpfr_init(r28091);
        mpfr_init(r28092);
        mpfr_init_set_str(r28093, "1", 10, MPFR_RNDN);
        mpfr_init(r28094);
        mpfr_init_set_str(r28095, "-4", 10, MPFR_RNDN);
        mpfr_init(r28096);
        mpfr_init(r28097);
        mpfr_init(r28098);
        mpfr_init(r28099);
        mpfr_init(r28100);
        mpfr_init(r28101);
        mpfr_init(r28102);
}

double f_fm(double wj, double x) {
        mpfr_set_d(r28090, wj, MPFR_RNDN);
        mpfr_mul(r28091, r28090, r28090, MPFR_RNDN);
        mpfr_sub(r28092, r28090, r28091, MPFR_RNDN);
        ;
        mpfr_div(r28094, r28093, r28090, MPFR_RNDN);
        ;
        mpfr_pow(r28096, r28094, r28095, MPFR_RNDN);
        mpfr_fma(r28097, r28090, r28092, r28096, MPFR_RNDN);
        mpfr_set_d(r28098, x, MPFR_RNDN);
        mpfr_exp(r28099, r28090, MPFR_RNDN);
        mpfr_fma(r28100, r28090, r28099, r28099, MPFR_RNDN);
        mpfr_div(r28101, r28098, r28100, MPFR_RNDN);
        mpfr_add(r28102, r28097, r28101, MPFR_RNDN);
        return mpfr_get_d(r28102, MPFR_RNDN);
}

static mpfr_t r28103, r28104, r28105, r28106, r28107, r28108, r28109, r28110, r28111, r28112, r28113, r28114, r28115;

void setup_mpfr_f_dm() {
        mpfr_set_default_prec(848);
        mpfr_init(r28103);
        mpfr_init(r28104);
        mpfr_init(r28105);
        mpfr_init_set_str(r28106, "1", 10, MPFR_RNDN);
        mpfr_init(r28107);
        mpfr_init_set_str(r28108, "-4", 10, MPFR_RNDN);
        mpfr_init(r28109);
        mpfr_init(r28110);
        mpfr_init(r28111);
        mpfr_init(r28112);
        mpfr_init(r28113);
        mpfr_init(r28114);
        mpfr_init(r28115);
}

double f_dm(double wj, double x) {
        mpfr_set_d(r28103, wj, MPFR_RNDN);
        mpfr_mul(r28104, r28103, r28103, MPFR_RNDN);
        mpfr_sub(r28105, r28103, r28104, MPFR_RNDN);
        ;
        mpfr_div(r28107, r28106, r28103, MPFR_RNDN);
        ;
        mpfr_pow(r28109, r28107, r28108, MPFR_RNDN);
        mpfr_fma(r28110, r28103, r28105, r28109, MPFR_RNDN);
        mpfr_set_d(r28111, x, MPFR_RNDN);
        mpfr_exp(r28112, r28103, MPFR_RNDN);
        mpfr_fma(r28113, r28103, r28112, r28112, MPFR_RNDN);
        mpfr_div(r28114, r28111, r28113, MPFR_RNDN);
        mpfr_add(r28115, r28110, r28114, MPFR_RNDN);
        return mpfr_get_d(r28115, MPFR_RNDN);
}