- iOS/Objective-C 元类和类别
- objective-c - -1001 错误,当 NSURLSession 通过 httpproxy 和/etc/hosts
- java - 使用网络类获取 url 地址
- ios - 推送通知中不播放声音
我正在使用 C++ 模拟涉及多个 (27) 刚性常微分方程的生物模型。我的程序在 MS C++ 2010 表达式编译器下完美运行,但在 g++ 编译器(NetBeans 6.8、Ubuntu 10.04 LTS)下运行失败。问题是一些变量变成了 NaN .以下是g++编译器下程序每一步后变量Vm
的值:
-59.4 -59.3993 -59.6081 100.081 34.6378 -50392.8 nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan nan
下面是相同代码在 MS C++ 编译器下未经任何更改的输出:
-59.4 -59.3993 -59.3986 -59.3979 -59.3972 -59.3966 -59.3959 -59.3952 -59.3946 -59.3939 -59.3933 -59.3926 -59.392 -59.3914 -59.3907 -59.3901 -59.3895 -59.3889 -59.3883 -59.3877 -59.3871 -59.3865 -59.3859 -59.3853 -59.3847 -59.3841 -59.3836 -59.383 -59.3824 -59.3819 -59.3813 -59.3808 -59.3802 -59.3797 -59.3791 -59.3786 -59.3781 -59.3775 -59.377 -
我只使用了“cmath”和“fstream”库。
问题出在哪里?两种场景的代码完全一样。
编辑 1:
好的,这是完整的代码:
#include <iostream>
#include <fstream>
#include <cmath>
using namespace std;
void FCN(void);
const int TMAX = 10000; //[ms] simulation time
const int TSTEP = 1;
const int MXSTEP = TMAX / TSTEP;
int ISTEPPRINT = MXSTEP / 5000; //time step for storing on disc
int ISTEP = 0;
const double R = 8341.0;
const double temp = 293.0;
const double F = 96487.0;
const double RT_F = R*temp / F;
const double z_K = 1;
const double z_Na = 1;
const double z_Ca = 2;
const double z_Cl = -1;
const double N_Av = 6.022e23;
double Ca_o = 2.0;
double Na_o = 140.0;
double Cl_o = 129.0;
double K_o = 5;
double NE = 0;
double NO = 0;
double cGMP = 0; //[mM] [cGMP]i
double cGMPprime = 0; //var
double IP3 = 0; //[mM] [IP3]i
double IP3d = 0; //var
double IP3prime = 0; //var
double DAG = 0; //[mM]
double DAGprime = 0; //var
double Ca_u = 0.66;
double Ca_r = 0.57;
double Ca_i = 68e-6;
double Na_i = 8.4;
double K_i = 140;
double Cl_i = 59.4;
double V_m = -59.4;
double V_mprime;
double Na_iprime;
double K_iprime;
double Cl_iprime;
double Ca_iprime;
double vol_i = 1;
double I_Natotm;
double I_Ktotm;
double I_Cltotm;
double I_Catotm; //[pA] total membrane Ca current
//Reversal potentials
double E_Ca; //[mV]
double E_Na; //[mV]
double E_K; //[mV]
double E_Cl;
//Membrane capacitance and area
double C_m = 25.0;
double A_m = C_m / 1e6;
//Voltage dependent calcium current I_CaL
double I_CaL;
double P_CaL = 1.88e-5;
double d_L;
double d_Lprime;
double d_Lbar;
double tau_d_L;
double f_L;
double f_Lbar;
double tau_f_L;
double f_Lprime;
//Delayed rectifier current I_K
double I_K;
double g_K = 1.35;
double p_K;
double p_Kbar;
double V_1_2 = -11.0;
double k = 15.0;
double tau_p_K;
double p_Kprime;
double q_1 = 1;
double q_2 = 1;
double q_bar;
double q_1prime;
double q_2prime;
double Pmin_NSC = 0.4344;
double Po_NSC;
double PNa_NSC = (5.11e-7);
double PCa_NSC = (5.11e-7)*4.54;
double PK_NSC = (5.11e-7)*1.06; //
double d_NSCmin = 0.0244;
double K_NSC = 3.0e-3;
double INa_NSC;
double ICa_NSC;
double IK_NSC;
double I_NSC;
//KATP current I_KATP
double I_KATP; //[pA] background K current
double g_KATP = 0.067; //[nS] max. background K current conductance
//Inward rectifier current I_K_i
double I_K_i; //[pA]
double g_maxK_i; //[nS] max. slope conductance
const double G_K_i = 0; // inward rectifier constant
const double n_K_i = 0.5; // inward rectifier constant
//Calcium-activated potassium current I_KCa
double I_KCa;
double i1_KCa;
double P_BKCa = 3.9e-13;
double N_BKCa = 6.6e6;
double P_KCa;
double p_obar;
double V_1_2_KCa;
double p_f;
double p_s;
double p_fprime;
double p_sprime;
double tau_pf = 0.84;
double tau_ps = 35.9;
double dV_1_2_KCa_NO = 46.3;
double R_NO;
double dV_1_2_KCa_cGMP = 76;
double R_cGMP;
double k_leak = 1;
double R_00;
double R_01 = 0.9955;
double R_10 = 0.0033;
double R_11 = 4.0e-6;
double R_01prime;
double R_10prime;
double R_11prime;
const double Kr1 = 2500.0;
const double Kr2 = 1.05;
const double K_r1 = 0.0076;
const double K_r2 = 0.084;
double I_up;
const double I_upbar = 3.34 * (k_leak + 1);
const double K_mup = 0.001;
double I_tr;
const double vol_u = 0.07;
double tau_tr = 1000.0;
double I_rel;
const double vol_r = 0.007;
const double tau_rel = 0.0333; //[ms]
const double R_leak = 1.07e-5 * (k_leak); ////equal to R_10^2 during concentration clamp
// time constant of SR release
double Ca_uprime; // dCa_u/dt
double Ca_rprime; // dCa_r/dt
double S_CM;
const double K_d = 2.6e-4;
const double S_CMbar = 0.1;
double CaCM;
const double K_dB = 5.298e-4;
const double B_Fbar = 0.1;
const double vol_Ca = 0.7;
const double CSQNbar = 15;
const double K_CSQN = 0.8;
double I_PMCA;
double I_PMCAbar = 5.37;
double K_mPMCA = 170e-6;
double I_NaK; ////[pA] Na/K pump
double I_NaKbar = 2.3083;
const double K_mK = 1.6;
const double K_mNa = 22;
const double Q_10_NaK = 1.87;
double I_NCX;
const double gamma2 = 0.45; //
double g_NCX = 0.000487; //[nS]
double d_NCX = 0.0003; //
double Fi_F; //
double Fi_R; //
double I_NaKCl_Cl; //[pA]
double L_cotr = 1.79e-8;
double I_M = 0; //[pA]
double I_MCa = 0;
double I_MNa = 0; //[pA] Na component
double I_MK = 0;
double I_SOC; //[pA]
double I_SOCCa;
double I_SOCNa; //[pA] Na component
const double g_SOCCa = 0.0083; //[nS]
const double g_SOCNa = 0.0575; //[nS]
const double H_SOC = 1;
const double K_SOC = 0.0001;
const double tau_SOC = 100;
double P_SOCbar;
double P_SOC = 0;
double P_SOCprime;
//Chloride currents
double I_Cl;
double g_Cl = 0.23;
double alpha_Cl;
double P_Cl;
//Stimulation current
double I_stim = 0; //[pA]
//IP3 receptor
double I_IP3;
double I_IP3bar = 2880e-6; //[1/ms]
double K_IP3 = 0.12e-3;
double K_actIP3 = 0.17e-3;
double K_inhIP3 = 0.1e-3; //[mM]
double h_IP3;
double k_onIP3 = 1.4;
double h_IP3prime;
double R_TG = 2e4;
double K_1G = 0.01;
double K_2G = 0.2;
double k_rG = 1.75e-7;
double k_pG = 0.1e-3;
double k_eG = 6e-6;
double ksi_G = 0.85;
double G_TG = 1e5;
double k_degG = 1.25e-3;
double k_aG = 0.17e-3;
double k_dG = 1.5e-3;
double PIP2_T = 5e7;
double r_rG = 0.015e-3;
double K_cG = 0.4e-3;
double alpha_G = 2.781e-8;
double vol_IP3 = vol_i;
double gamma_G = N_Av*vol_IP3 * 1e-15;
double RS_G = R_TG*ksi_G;
double RS_PG = 0;
double PIP2; //
double r_hG;
double G;
double delta_G; //
double RS_Gprime;
double RS_PGprime;
double Gprime;
double PIP2prime;
double rho_rG;
//cGMP formation
double k1sGC = 2e3; //[1/mM/ms]
double k_1sGC = 15e-3; //[1/ms]
double k2sGC = 0.64e-5; //[1/ms]
double k_2sGC = 0.1e-6; //[1/ms]
double k3sGC = 4.2; //[1/mM/ms]
double kDsGC = 0.4e-3;
double kDact_deactsGC = 0.1e-3; //[1/ms]
double V_cGMP = 0; //
double V_cGMPprime;
double V_cGMPmax = 0.1 * 1.26e-6; //[mM/ms]
double V_cGMPbar;
double B5sGC = k2sGC / k3sGC;
double A0sGC = ((k_1sGC + k2sGC) * kDsGC + k_1sGC*k_2sGC) / (k1sGC*k3sGC);
double A1sGC = ((k1sGC + k3sGC) * kDsGC + (k2sGC + k_2sGC) * k1sGC) / (k1sGC*k3sGC);
double kpde_cGMP = 0.0695e-3; //[1/ms]
double tausGC;
const int N = 27;
double Y[N], Y1[N], YPRIME[N];
int N1 = 33;
double T = 0;
int main(void) {
ofstream fileY, fileY1, fileT;
// initial conditions SMC
//ICaL
d_Lbar = 1.0 / (1 + exp(-(V_m) / 8.3));
d_L = d_Lbar;
f_Lbar = 1.0 / (1 + exp((V_m + 42.0) / 9.1));
f_L = f_Lbar;
//IKCa
R_NO = NO / (NO + 200e-6);
R_cGMP = pow(cGMP, 2) / (pow(cGMP, 2) + pow(0.55e-3, 2));
V_1_2_KCa = -41.7 * log10(Ca_i) - 128.2 - dV_1_2_KCa_NO * R_NO - dV_1_2_KCa_cGMP*R_cGMP;
p_obar = 1 / (1 + exp(-(V_m - V_1_2_KCa) / 18.25));
p_f = p_obar;
p_s = p_obar;
//I_K
p_Kbar = 1 / (1 + exp(-(V_m - V_1_2) / k));
p_K = p_Kbar;
q_bar = 1.0 / (1 + exp((V_m + 40) / 14));
q_1 = q_bar;
q_2 = q_bar;
//IP3 receptor
h_IP3 = K_inhIP3 / (Ca_i + K_inhIP3);
//norepinephrine receptor
PIP2 = PIP2_T - (1 + k_degG / r_rG) * gamma_G*IP3;
r_hG = k_degG * gamma_G * IP3 / PIP2;
G = (K_cG + Ca_i) / (alpha_G * Ca_i) * r_hG;
delta_G = k_dG * G / (k_aG * (G_TG - G));
Y[0] = V_m;
Y[1] = d_L;
Y[2] = f_L;
Y[3] = p_K;
Y[4] = q_1;
Y[5] = p_f;
Y[6] = p_s;
Y[7] = R_01;
Y[8] = R_10;
Y[9] = R_11;
Y[10] = Ca_u;
Y[11] = Ca_r;
Y[12] = Ca_i;
Y[13] = Na_i;
Y[14] = K_i;
Y[15] = q_2;
Y[16] = P_SOC;
Y[17] = Cl_i;
Y[18] = h_IP3;
Y[19] = RS_G;
Y[20] = RS_PG;
Y[21] = G;
Y[22] = IP3;
Y[23] = PIP2;
Y[24] = DAG;
Y[25] = V_cGMP;
Y[26] = cGMP;
ISTEP = -1;
T = 0.0 - TSTEP;
fileY.open("Y.txt");
fileY1.open("Y1.txt");
fileT.open("T.txt");
for (;;) {
ISTEP = ISTEP + 1;
T = T + TSTEP;
//Norepinephrine
if (T > 10000) NE = 1e-3; //NE [mM] beginning of stimulation
if (T > 70000) NE = 0; //end of stimulation
//Nitric oxide
//IF (T>30000) NO = 1D-3 //NO [mM]
//IF (T>70000) NO = 0
//Extracellular potassium
//IF (T>10000) K_o = 30
//IF (T>70000) K_o = 5
//Current
//IF (T>10000) I_stim = 5 //I_stim [pA] current injection
//IF (T>40000) I_stim = -5
//IF (T>70000) I_stim = 0
// For the time being, I just interested in Y[0] values (which is V_m actually)
fileY << Y[0];
fileY << "\t";
if ((ISTEP % ISTEPPRINT) == 0) {
// for (int i=0; i< N; i++) {
// fileY << Y[i];
// fileY << "\t";
// }
// fileY << endl;
// for (int i=0; i< N1; i++) {
// fileY1 << Y1[i];
// fileY1 << "\t";
// }
// fileY1 << endl;
//
//
//
// fileT << T;
// fileT << "\t";
}
FCN();
for (int i = 0; i < N; i++) {
Y[i] = Y[i] + TSTEP * YPRIME[i];
}
// disp(Yconcat(1))
if (ISTEP == MXSTEP)
break;
}
cout << "It is done!" << endl;
cout << Y[0] << endl;
fileY.close();
fileY1.close();
fileT.close();
return 0;
}
void FCN(void) {
V_m = Y[0];
d_L = Y[1];
f_L = Y[2];
p_K = Y[3];
q_1 = Y[4];
p_f = Y[5];
p_s = Y[6];
R_01 = Y[7];
R_10 = Y[8];
R_11 = Y[9];
Ca_u = Y[10];
Ca_r = Y[11];
Ca_i = Y[12];
Na_i = Y[13];
K_i = Y[14];
q_2 = Y[15];
P_SOC = Y[16];
Cl_i = Y[17];
h_IP3 = Y[18];
RS_G = Y[19];
RS_PG = Y[20];
G = Y[21];
IP3 = Y[22];
PIP2 = Y[23];
DAG = Y[24];
V_cGMP = Y[25];
cGMP = Y[26];
//-------------------------------------- Model equations ---------------------------------------------
//cGMP formation
V_cGMPbar = V_cGMPmax * (B5sGC * NO + pow(NO, 2)) / (A0sGC + A1sGC * NO + pow(NO, 2));
if ((V_cGMPbar - V_cGMP) >= 0) {
tausGC = 1 / (k3sGC * NO + kDact_deactsGC); //kDact_deactsGC different from original kDsGC to uncouple Km from time constants
} else {
tausGC = 1 / (kDact_deactsGC + k_2sGC);
}
V_cGMPprime = (V_cGMPbar - V_cGMP) / tausGC;
cGMPprime = V_cGMP - kpde_cGMP * cGMP * cGMP / (1e-3 + cGMP);
//Norepinephrine receptor
RS_Gprime = (k_rG * ksi_G * R_TG - (k_rG + k_pG * NE / (K_1G + NE)) * RS_G - k_rG * RS_PG);
RS_PGprime = NE * (k_pG * RS_G / (K_1G + NE) - k_eG * RS_PG / (K_2G + NE));
rho_rG = NE * RS_G / (ksi_G * R_TG * (K_1G + NE));
Gprime = k_aG * (delta_G + rho_rG)*(G_TG - G) - k_dG*G;
r_hG = alpha_G * Ca_i / (K_cG + Ca_i) * G;
IP3prime = r_hG / gamma_G * PIP2 - k_degG*IP3;
PIP2prime = -(r_hG + r_rG) * PIP2 - r_rG * gamma_G * IP3 + r_rG*PIP2_T;
DAGprime = r_hG / gamma_G * PIP2 - k_degG*DAG;
//Reversal potentials
E_Ca = RT_F / z_Ca * log(Ca_o / Ca_i);
E_Na = RT_F * log(Na_o / Na_i);
E_K = RT_F * log(K_o / K_i);
E_Cl = RT_F / z_Cl * log(Cl_o / Cl_i);
//Voltage dependent calcium current I_CaL
tau_d_L = 2.5 * exp(-pow((V_m + 40) / 30, 2)) + 1.15;
d_Lbar = 1.0 / (1 + exp(-(V_m) / 8.3));
d_Lprime = (d_Lbar - d_L) / tau_d_L;
f_Lbar = 1.0 / (1 + exp((V_m + 42.0) / 9.1));
tau_f_L = 65 * exp(-pow((V_m + 35) / 25, 2)) + 45;
f_Lprime = (f_Lbar - f_L) / tau_f_L;
if (V_m == 0) {
I_CaL = d_L * f_L * P_CaL * A_m * 1e6 * z_Ca * F * (Ca_i - Ca_o); //[pA]
} else {
I_CaL = d_L * f_L * P_CaL * A_m * 1e6 * V_m * pow(z_Ca * F, 2) / (R * temp)*(Ca_o - Ca_i * exp(V_m * z_Ca / (RT_F))) / (1 - exp(V_m * z_Ca / (RT_F))); //[pA]
}
//Delayed rectifier current I_K
p_Kbar = 1 / (1 + exp(-(V_m - V_1_2) / k));
tau_p_K = 61.49 * exp(-0.0268 * V_m);
p_Kprime = (p_Kbar - p_K) / tau_p_K;
q_bar = 1.0 / (1 + exp((V_m + 40) / 14));
q_1prime = (q_bar - q_1) / 371;
q_2prime = (q_bar - q_2) / 2884;
I_K = 1 * g_K * p_K * (0.45 * q_1 + 0.55 * q_2) * (V_m - E_K);
//Alpha-adrenoceptor-activated nonselective cation channel NSC
Po_NSC = Pmin_NSC + (1 - Pmin_NSC) / (1 + exp(-(V_m - 47.12) / 24.24));
if (V_m == 0) {
INa_NSC = 1 * (DAG / (DAG + K_NSC) + d_NSCmin) * Po_NSC * PNa_NSC * A_m * 1e6 * F * (Na_i - Na_o);
ICa_NSC = 1 * (0 * DAG / (DAG + K_NSC) + d_NSCmin) * Po_NSC * PCa_NSC * A_m * 1e6 * z_Ca * F * (Ca_i - Ca_o);
IK_NSC = 1 * (DAG / (DAG + K_NSC) + d_NSCmin) * Po_NSC * PK_NSC * A_m * 1e6 * F * (K_i - K_o);
} else {
INa_NSC = 1 * (DAG / (DAG + K_NSC) + d_NSCmin) * Po_NSC * PNa_NSC * A_m * 1e6 * V_m * pow(F, 2) / (R * temp)*(Na_o - Na_i * exp(V_m / RT_F)) / (1 - exp(V_m / RT_F));
ICa_NSC = 1 * (0 * DAG / (DAG + K_NSC) + d_NSCmin) * Po_NSC * PCa_NSC * A_m * 1e6 * V_m * pow(z_Ca * F, 2) / (R * temp)*(Ca_o - Ca_i * exp(V_m * z_Ca / RT_F)) / (1 - exp(V_m * z_Ca / RT_F));
IK_NSC = 1 * (DAG / (DAG + K_NSC) + d_NSCmin) * Po_NSC * PK_NSC * A_m * 1e6 * V_m * pow(F, 2) / (R * temp)*(K_o - K_i * exp(V_m / RT_F)) / (1 - exp(V_m / RT_F));
}
I_NSC = ICa_NSC + INa_NSC + IK_NSC;
//KATP current I_KATP
I_KATP = g_KATP * (V_m - E_K);
//Inward rectifier current I_K_i
g_maxK_i = G_K_i * pow(K_o, n_K_i);
I_K_i = g_maxK_i * (V_m - E_K) / (1 + exp((V_m - E_K) / 28.89));
//Calcium-activated potassium current I_KCa
if (V_m == 0) {
i1_KCa = 1e6 * P_BKCa * F * (K_i - K_o); //[pA]
} else {
i1_KCa = 1e6 * P_BKCa * V_m * F / RT_F * (K_o - K_i * exp(V_m / RT_F)) / (1 - exp(V_m / RT_F)); //[pA]
} //Mistry and Garland 1998
R_NO = NO / (NO + 200e-6);
R_cGMP = pow(cGMP, 2) / (pow(cGMP, 2) + pow(1.5e-3, 2));
V_1_2_KCa = -41.7 * log10(Ca_i) - 128.2 - dV_1_2_KCa_NO * R_NO - dV_1_2_KCa_cGMP*R_cGMP;
p_obar = 1 / (1 + exp(-(V_m - V_1_2_KCa) / 18.25));
p_fprime = (p_obar - p_f) / tau_pf;
p_sprime = (p_obar - p_s) / tau_ps;
P_KCa = 0.17 * p_f + 0.83 * p_s;
I_KCa = A_m * N_BKCa * i1_KCa * P_KCa;
//Store operated non-selective cation channel
P_SOCbar = 1 / (1 + pow(Ca_u / K_SOC, H_SOC));
P_SOCprime = (P_SOCbar - P_SOC) / tau_SOC;
I_SOCCa = 1 * P_SOC * g_SOCCa * (V_m - E_Ca);
I_SOCNa = 1 * P_SOC * g_SOCNa * (V_m - E_Na);
I_SOC = I_SOCCa + I_SOCNa;
//Chloride currents
alpha_Cl = pow(cGMP, 3.3) / (pow(cGMP, 3.3) + pow(6.4e-3, 3.3));
P_Cl = pow(Ca_i, 2) / (pow(Ca_i, 2) + pow(365e-6, 2)) * 0.0132 + pow(Ca_i, 2) / (pow(Ca_i, 2) + pow(400e-6 * (1 - alpha_Cl * 0.9), 2)) * alpha_Cl;
I_Cl = P_Cl * g_Cl * C_m * (V_m - E_Cl);
//IP3 receptor current
h_IP3prime = k_onIP3 * (K_inhIP3 - (Ca_i + K_inhIP3) * h_IP3);
I_IP3 = I_IP3bar * pow(IP3 / (IP3 + K_IP3) * Ca_i / (Ca_i + K_actIP3) * h_IP3, 3)*(Ca_u - Ca_i) * z_Ca * F*vol_Ca;
//Calcium-induced calcium release
R_00 = 1 - R_01 - R_10 - R_11;
R_10prime = Kr1 * pow(Ca_i, 2) * R_00 - (K_r1 + Kr2 * Ca_i) * R_10 + K_r2*R_11;
R_11prime = Kr2 * Ca_i * R_10 - (K_r1 + K_r2) * R_11 + Kr1 * pow(Ca_i, 2) * R_01;
R_01prime = Kr2 * Ca_i * R_00 + K_r1 * R_11 - (K_r2 + Kr1 * pow(Ca_i, 2)) * R_01;
I_up = I_upbar * Ca_i / (Ca_i + K_mup);
I_tr = (Ca_u - Ca_r) * (2 * F * vol_u) / tau_tr;
I_rel = (pow(R_10, 2) + R_leak) * (Ca_r - Ca_i) * (2 * F * vol_r) / tau_rel;
Ca_uprime = (I_up - I_tr - I_IP3) / (2 * F * vol_u);
Ca_rprime = (I_tr - I_rel) / (2 * F * vol_r) / (1 + CSQNbar * K_CSQN / pow((K_CSQN + Ca_r), 2));
//Ca buffering and cytosolic material balance
S_CM = S_CMbar * K_d / (K_d + Ca_i);
CaCM = S_CMbar - S_CM;
I_PMCA = I_PMCAbar * Ca_i / (Ca_i + K_mPMCA);
//NaK pump
I_NaK = pow(Q_10_NaK, ((temp - 309.15) / 10)) * C_m * I_NaKbar * ((pow(K_o, 1.1)) / (pow(K_o, 1.1) + (pow(K_mK, 1.1)))
*(pow(Na_i, 1.7)) / ((pow(Na_i, 1.7))+(pow(K_mNa, 1.7)))) *(V_m + 150) / (V_m + 200);
Fi_F = exp(gamma2 * V_m * F / (R * temp));
Fi_R = exp((gamma2 - 1) * V_m * F / (R * temp));
I_NCX = 1 * (1 + 0.55 * cGMP / (cGMP + (45e-3))) * g_NCX * (pow(Na_i, 3) * Ca_o * Fi_F - pow(Na_o, 3) * Ca_i * Fi_R) / (1 + d_NCX * (pow(Na_o, 3) * Ca_i + pow(Na_i, 3) * Ca_o));
I_NaKCl_Cl = (1 + 7 / 2 * cGMP / (cGMP + 6.4e-3))*(-A_m * L_cotr * R * temp * z_Cl * F * log(Na_o / Na_i * K_o / K_i * pow(Cl_o / Cl_i, 2)));
I_Catotm = I_SOCCa + I_CaL - 2 * I_NCX + I_PMCA + ICa_NSC + I_MCa;
Ca_iprime = -(I_Catotm + I_up - I_rel - I_IP3) / (2 * F * vol_Ca) / (1 + S_CMbar * K_d / (pow(K_d + Ca_i, 2)) + B_Fbar * K_dB / (pow(K_dB + Ca_i, 2)));
I_Natotm = -0.5 * I_NaKCl_Cl + I_SOCNa + 3 * I_NaK + 3 * I_NCX + INa_NSC + I_MNa;
Na_iprime = -(I_Natotm) / (F * vol_i);
I_Ktotm = -0.5 * I_NaKCl_Cl + I_K + I_KCa + I_K_i + IK_NSC + I_KATP - 2 * I_NaK + I_MK;
K_iprime = -(I_Ktotm) / (F * vol_i);
I_Cltotm = I_NaKCl_Cl + I_Cl;
Cl_iprime = -(I_Cltotm) / (z_Cl * F * vol_i);
//Transmembrane potential
V_mprime = -1 / C_m * (-I_stim + I_Cl + I_SOC + I_CaL + I_K + I_KCa + I_K_i + I_M + I_NCX + I_NaK + I_PMCA + I_NSC + I_KATP);
//YPRIME = zeros(1, N);
YPRIME[0] = V_mprime;
YPRIME[1] = d_Lprime;
YPRIME[2] = f_Lprime;
YPRIME[3] = p_Kprime;
YPRIME[4] = q_1prime;
YPRIME[5] = p_fprime;
YPRIME[6] = p_sprime;
YPRIME[7] = R_01prime;
YPRIME[8] = R_10prime;
YPRIME[9] = R_11prime;
YPRIME[10] = Ca_uprime;
YPRIME[11] = Ca_rprime;
YPRIME[12] = Ca_iprime;
YPRIME[13] = Na_iprime;
YPRIME[14] = K_iprime;
YPRIME[15] = q_2prime;
YPRIME[16] = P_SOCprime;
YPRIME[17] = Cl_iprime;
YPRIME[18] = h_IP3prime;
YPRIME[19] = RS_Gprime;
YPRIME[20] = RS_PGprime;
YPRIME[21] = Gprime;
YPRIME[22] = IP3prime;
YPRIME[23] = PIP2prime;
YPRIME[24] = DAGprime;
YPRIME[25] = V_cGMPprime;
YPRIME[26] = cGMPprime;
//Non state variables
Y1[0] = I_CaL;
Y1[1] = I_K;
Y1[2] = I_K_i;
Y1[3] = I_NSC;
Y1[4] = I_KCa;
Y1[5] = I_up;
Y1[6] = I_rel;
Y1[7] = I_PMCA;
Y1[8] = I_NCX;
Y1[9] = I_NaK;
Y1[10] = INa_NSC;
Y1[11] = ICa_NSC;
Y1[12] = IK_NSC;
Y1[13] = I_SOCCa;
Y1[14] = I_SOCNa;
Y1[15] = I_Cl;
Y1[16] = I_NaKCl_Cl;
Y1[17] = I_IP3;
Y1[18] = I_tr;
Y1[19] = NE;
Y1[20] = I_KATP;
Y1[21] = I_stim;
Y1[22] = V_cGMPbar;
Y1[23] = NO;
Y1[29] = I_Catotm;
Y1[30] = I_Natotm;
Y1[31] = I_Cltotm;
Y1[32] = I_Ktotm;
}
最佳答案
我怀疑您的代码正在为您的 MS C 编译器提供答案,但我怀疑它能否完美运行。 NaN(非数字)是计算超出其范围的函数的结果。因为你没有提到你是如何解决你的僵硬系统的,我不知道你是在取负数的对数,还是其他一些松鼠算术,但我确信你在你的 Windows 代码中做同样的算术。仅仅因为代码试图蒙混过关并不意味着它正在正确运行。
我建议您开始查看这两个程序从哪里开始出现分歧,您可能会发现一个有问题的操作,即 MS 编译器在 g++ 返回 NaN 的情况下返回 0 或 Inf。
关于c++ - 我的 C++ 代码在 MS C++ 编译器上运行完美,但在 g++ 编译器上给我 NaN。为什么?,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/3295856/
我已经使用 MS Word 和一大堆表单字段创建了一个应用程序表单,并且我有一个 Access db,可以从这个 Word 文档中导入我需要的所有数据,这要归功于: http://msdn.micro
我试图找到一种将 Outlook 插件发布到办公商店的方法。但我发现我们只能发布 Office 应用程序,而不能发布 Office 商店的加载项。因此我想知道 Office 应用程序和 Office
我在 MS Reporting Services 服务器上部署了一份报告,工作正常。我可以使用 Microsoft 的报表查看器组件从 ASPX 页面毫无问题地 Access 它、设置报表参数等。效果
让我们再试一次。我发布这个是为了回答 2 个问题 MS Project 2007 是否需要 SharePoint(我希望没有)? 做 你喜欢 MS Project 开发团队 - 它是有用的还是 疼痛?
我正在执行这些星期六上午的任务之一,试图理解为什么为什么要在计算机注册表中搜索某些信息会花费大量时间,甚至迫使我停止该过程。使用这些注册表清理程序之一,我发现该代码花了数十分钟遍历如下行: HKEY_
从多年前开始,我就没有使用Access。 它能很好地解决什么样的问题,甚至比真正的RDBMS支持的Web应用程序更好? 它仍在积极开发吗?还是MS已经死了? 最大的局限性是什么? 更新: 应该使用什么
我们计划重新设计一个相当庞大的 MS Access 应用程序。有没有办法在同一应用程序上同时工作,或者是否可以合并同一文件的两个单独实例(不是数据,而是表单和代码)。现在 Access 包含数据,但在
我写了一些SQL命令来更正表中的字段。由于它是如此之小(也许我有点自大),我什至没有运行过一次,只是将其放入了更新包中供其他用户使用。 Dim SQL As String Dim rs As DAO.
它是Office自带的,是一个“中规中矩”的数据库,到今天这里有800多个问题,但我从来没有关注过它。 我失去了一些有趣的东西? 我说的是 MS-Access 作为用于快速原型(prototype)制
我有一个MS-Access数据库,该数据库已通过使用“用于Access的Microsoft SQL Server迁移助手2008”(aka SSMA)转换为使用SQL表并创建了链接表(因此,MS-Ac
我有一个 Excel 文件,其中包含从 Access 数据库(主数据库)导出的任务。然后,此 Excel 文件用作 MS Project 的导入文件。随后,MS Project 用于实际跟踪和报告,并
我正在尝试获取有关如何将 MS Project 2010 连接到 MS Project Server 2010 的教程或分步说明。 我已经在我的服务器上安装了 Server 2008 R2(64 位)
有没有办法像选择查询一样在查询中引用表单的组合框/文本框? 我通常在选择查询的条件中使用类似这样的东西: like forms!frmMain.qTitleofSomething&* (acces
我想创建一个表,其中包含 DOUBLE 实数类型的列。我可以在表设计 View 中找到数据类型 Number,但是没有 Double 或 single,Float.. 如何实现..? 我还需要 SQL
我环顾四周,发现了一些关于如何从字段的“描述”框中获取描述的 VBA 代码,但没有找到如何在表单属性中使用它的方法。 我希望出现一个 ControlTip,其中包含从数据库中的描述中带来的该字段的描述
我有一个难题。我已经开发了一个 Access 应用程序,我正准备分发它。我刚刚拆分了数据库。 (我知道,有人说我应该从一开始就把它分开开发……我没有)我也刚刚加密了后端数据库。在前端,我已链接到后端并
我制作了一个 MS Access 2013 数据库来跟踪有关交易网站的所有通信。与此问题相关的表和列是具有列 ID(编号)和链接(超链接)的广告,以及具有列广告的注释,其中包含广告 ID。链接字段包含
我与我不拥有且无法更改的数据库建立了 ODBC 连接。我要做的是使相关记录合并为一条记录。关系是一对多。 我有一个学生管理系统,想要导出一个提供自动标注服务(由调用收费)的调用列表。如果有多个学生住在
我在 Access 的表单中设置了一个文本框。该表单链接到一个表格。但是它自己的文本框是未绑定(bind)的,它用于简单地收集用户输入。但是,我无法编辑它所查看的值。 文本框未锁定。文本框可以在 VB
很难说出这里问的是什么。这个问题是模棱两可的、模糊的、不完整的、过于宽泛的或修辞的,无法以目前的形式得到合理的回答。如需帮助澄清这个问题以便重新打开它,visit the help center .
我是一名优秀的程序员,十分优秀!