double f(double c_p, double c_n, double t, double s) {
double r8519705 = 1.0;
double r8519706 = s;
double r8519707 = -r8519706;
double r8519708 = exp(r8519707);
double r8519709 = r8519705 + r8519708;
double r8519710 = r8519705 / r8519709;
double r8519711 = c_p;
double r8519712 = pow(r8519710, r8519711);
double r8519713 = r8519705 - r8519710;
double r8519714 = c_n;
double r8519715 = pow(r8519713, r8519714);
double r8519716 = r8519712 * r8519715;
double r8519717 = t;
double r8519718 = -r8519717;
double r8519719 = exp(r8519718);
double r8519720 = r8519705 + r8519719;
double r8519721 = r8519705 / r8519720;
double r8519722 = pow(r8519721, r8519711);
double r8519723 = r8519705 - r8519721;
double r8519724 = pow(r8519723, r8519714);
double r8519725 = r8519722 * r8519724;
double r8519726 = r8519716 / r8519725;
return r8519726;
}