#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 r20043 = wj;
        float r20044 = exp(r20043);
        float r20045 = r20043 * r20044;
        float r20046 = x;
        float r20047 = r20045 - r20046;
        float r20048 = r20044 + r20045;
        float r20049 = r20047 / r20048;
        float r20050 = r20043 - r20049;
        return r20050;
}

double f_id(double wj, double x) {
        double r20051 = wj;
        double r20052 = exp(r20051);
        double r20053 = r20051 * r20052;
        double r20054 = x;
        double r20055 = r20053 - r20054;
        double r20056 = r20052 + r20053;
        double r20057 = r20055 / r20056;
        double r20058 = r20051 - r20057;
        return r20058;
}


double f_of(float wj, float x) {
        float r20059 = wj;
        float r20060 = r20059 * r20059;
        float r20061 = r20059 * (r20059 * r20059);
        float r20062 = r20060 - r20061;
        float r20063 = x;
        float r20064 = 1.0f;
        float r20065 = r20059 + r20064;
        float r20066 = exp(r20059);
        float r20067 = r20065 * r20066;
        float r20068 = r20063 / r20067;
        float r20069 = r20062 + r20068;
        return r20069;
}

double f_od(double wj, double x) {
        double r20070 = wj;
        double r20071 = r20070 * r20070;
        double r20072 = r20070 * (r20070 * r20070);
        double r20073 = r20071 - r20072;
        double r20074 = x;
        double r20075 = 1.0;
        double r20076 = r20070 + r20075;
        double r20077 = exp(r20070);
        double r20078 = r20076 * r20077;
        double r20079 = r20074 / r20078;
        double r20080 = r20073 + r20079;
        return r20080;
}

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 r20081, r20082, r20083, r20084, r20085, r20086, r20087, r20088;

void setup_mpfr_f_im() {
        mpfr_set_default_prec(144);
        mpfr_init(r20081);
        mpfr_init(r20082);
        mpfr_init(r20083);
        mpfr_init(r20084);
        mpfr_init(r20085);
        mpfr_init(r20086);
        mpfr_init(r20087);
        mpfr_init(r20088);
}

double f_im(double wj, double x) {
        mpfr_set_d(r20081, wj, MPFR_RNDN);
        mpfr_exp(r20082, r20081, MPFR_RNDN);
        mpfr_mul(r20083, r20081, r20082, MPFR_RNDN);
        mpfr_set_d(r20084, x, MPFR_RNDN);
        mpfr_sub(r20085, r20083, r20084, MPFR_RNDN);
        mpfr_add(r20086, r20082, r20083, MPFR_RNDN);
        mpfr_div(r20087, r20085, r20086, MPFR_RNDN);
        mpfr_sub(r20088, r20081, r20087, MPFR_RNDN);
        return mpfr_get_d(r20088, MPFR_RNDN);
}

static mpfr_t r20089, r20090, r20091, r20092, r20093, r20094, r20095, r20096, r20097, r20098, r20099;

void setup_mpfr_f_fm() {
        mpfr_set_default_prec(144);
        mpfr_init(r20089);
        mpfr_init(r20090);
        mpfr_init(r20091);
        mpfr_init(r20092);
        mpfr_init(r20093);
        mpfr_init_set_str(r20094, "1", 10, MPFR_RNDN);
        mpfr_init(r20095);
        mpfr_init(r20096);
        mpfr_init(r20097);
        mpfr_init(r20098);
        mpfr_init(r20099);
}

double f_fm(double wj, double x) {
        mpfr_set_d(r20089, wj, MPFR_RNDN);
        mpfr_sqr(r20090, r20089, MPFR_RNDN);
        mpfr_mul(r20091, r20089, r20089, MPFR_RNDN); mpfr_mul(r20091, r20091, r20089, MPFR_RNDN);
        mpfr_sub(r20092, r20090, r20091, MPFR_RNDN);
        mpfr_set_d(r20093, x, MPFR_RNDN);
        ;
        mpfr_add(r20095, r20089, r20094, MPFR_RNDN);
        mpfr_exp(r20096, r20089, MPFR_RNDN);
        mpfr_mul(r20097, r20095, r20096, MPFR_RNDN);
        mpfr_div(r20098, r20093, r20097, MPFR_RNDN);
        mpfr_add(r20099, r20092, r20098, MPFR_RNDN);
        return mpfr_get_d(r20099, MPFR_RNDN);
}

static mpfr_t r20100, r20101, r20102, r20103, r20104, r20105, r20106, r20107, r20108, r20109, r20110;

void setup_mpfr_f_dm() {
        mpfr_set_default_prec(144);
        mpfr_init(r20100);
        mpfr_init(r20101);
        mpfr_init(r20102);
        mpfr_init(r20103);
        mpfr_init(r20104);
        mpfr_init_set_str(r20105, "1", 10, MPFR_RNDN);
        mpfr_init(r20106);
        mpfr_init(r20107);
        mpfr_init(r20108);
        mpfr_init(r20109);
        mpfr_init(r20110);
}

double f_dm(double wj, double x) {
        mpfr_set_d(r20100, wj, MPFR_RNDN);
        mpfr_sqr(r20101, r20100, MPFR_RNDN);
        mpfr_mul(r20102, r20100, r20100, MPFR_RNDN); mpfr_mul(r20102, r20102, r20100, MPFR_RNDN);
        mpfr_sub(r20103, r20101, r20102, MPFR_RNDN);
        mpfr_set_d(r20104, x, MPFR_RNDN);
        ;
        mpfr_add(r20106, r20100, r20105, MPFR_RNDN);
        mpfr_exp(r20107, r20100, MPFR_RNDN);
        mpfr_mul(r20108, r20106, r20107, MPFR_RNDN);
        mpfr_div(r20109, r20104, r20108, MPFR_RNDN);
        mpfr_add(r20110, r20103, r20109, MPFR_RNDN);
        return mpfr_get_d(r20110, MPFR_RNDN);
}

