Hi,
I am solving a large nonlinear system of equations. It contains 44 variables and equations.
But I face the error :
ExtraVariablesSystemException: The system is unbalanced. There are 47 highest order derivative variables and 44 equations.
using Pkg
Pkg.add(["DataFrames", "XLSX"])
Pkg.add("Optim")
##############Reading data from excel #################
using DataFrames, XLSX
df = DataFrame(XLSX.readtable("??????\SolverJ\\psc_ft_scenario_new.xlsx", 1))
i = 1 # Julia indexing starts at 1, reading te first row
Qsl = df[i, "Qsl"]
Qdil = df[i, "Qdil"]
Qmo = df[i, "Qm"]
Qto = df[i, "Qt"]
Qft = df[i, "Qft"]
α_bsl_1 = df[i, "a_bsl_1"]
α_bsl_2 = df[i, "a_bsl_2"]
α_bsl_3 = df[i, "a_bsl_3"]
α_bsl_4 = df[i, "a_bsl_4"]
α_ssl_1 = df[i, "a_ssl_1"]
α_ssl_2 = df[i, "a_ssl_2"]
α_ssl_3 = df[i, "a_ssl_3"]
Qfd = Qsl + Qdil
# PSV FEED COMPOSITION
α_bfd_1 = α_bsl_1 * Qsl / Qfd
α_bfd_2 = α_bsl_2 * Qsl / Qfd
α_bfd_3 = α_bsl_3 * Qsl / Qfd
α_bfd_4 = α_bsl_4 * Qsl / Qfd
α_sfd_1 = α_ssl_1 * Qsl / Qfd
α_sfd_2 = α_ssl_2 * Qsl / Qfd
α_sfd_3 = α_ssl_3 * Qsl / Qfd
α_wfd = 1 - (α_bfd_1 + α_bfd_2 + α_bfd_3 + α_bfd_4 + α_sfd_1 + α_sfd_2 + α_sfd_3)
# ITERATIVELY RETRIEVE THE DATA FOR DIFFERENT SCENARIO
###########################################################
using NLsolve, NonlinearSolve, Plots, ModelingToolkit
@variables begin
v_bf_1 = -0.04660483043350957
v_bf_2 = -0.09345981917569582
v_sf_1 = -0.08130625749545124
v_bmu_1 = -0.09613872404222351
v_bmu_2 = -0.030677683191747712
v_smu_1 = -0.07997942776010468
v_bt_1 = 0.013060048582814041
v_st_1 = 0.05644695757124584
v_st_2 = 0.09751454192570927
v_st_3 = 0.08123282460225396
v_wf = -0.058046308842861466
v_wmu = -0.00036375845116324483
v_wt = 0.007696945147263263
#####
ρ_f = 960.8014204117135
ρ_mu = 975.0320325950338
ρ_m = 1138.387734137675
ρ_t= 1476.4611449382703
#######
α_bf = 0.7397525076061395 # 0.6397525076
α_sf = 0.05714768930793848 # 0.57147689307938
α_bmu = 0.6852046880649323
α_smu = 0.053434669298424966 # 0.43434669298424966
α_bm = 0.07531046509265504 # makes sense
α_sm = 0.21740611040100571
α_bt = 0.09304537667828224 # 0.49008855975103
α_st = 0.49008855975103 # 0.09304537667828224
α_bf_1 = 0.018487075313053238 # 0.37705165301787147
α_bf_2 = 0.711263232293 # 0.018487075313053238
α_sf_1 = 0.040193096
α_bmu_1 = 0.2058157585672766 # 0.4758001863282879
α_bmu_2 = 0.4758001863282879 # 0.9258157585672766
α_smu_1 = 0.05702341 # 0.8259463792344883
α_bm_1 = 0.010464741819863431 # 0.10464741819863431
α_bm_2 = 0.06144341881817391 # 0.08144341881817391
α_sm_1 = 0.10082698 # 0.14524751024977683
α_sm_2 = 0.010710206940301835 # 0.1871020694030183
α_sm_3 = 0.10586891754887175 # 0.2586891754887175
α_bt_1 = 0.0045 # 0.11805236842596298
α_st_1 = 0.01384434783970221
α_st_2 = 0.042886064930470547 # 0.42886064930470547
α_st_3 = 0.60029809495899 # 0.6190129809495899
###
α_wf = 0.20309981069 # 0.4277137177198411
α_wmu = 0.2613606426 # 0.9023711547833565
α_wm = 0.7072834245 # changed with the previous value and its interchanged
α_wt = 0.3384714 # 0.4533852950765186
###
VI = -9.541152104376794e-11
###
w_bf = 0.6403761770580065 # 0.8403761770580065
w_sf = 0.1299167398195812 # 0.3799167398195812
w_wf = 0.255356814994774 # 0.735356814994774
w_bm = 0.031998081785725 # 0.831998081785725
w_sm = 0.1946105277644868 # 0.7946105277644868
w_wm = 0.7622801692643188 # 0.5722801692643188
w_bt = 0.0442609239538598 # 0.44426092395385985
w_st = 0.5819373884863865 # 0.9819373884863865
w_wt = 0.2737807493701219 # 0.2737807493701219
###
Qff = 0.399794587515803
Qf = 0.5301731026307072
Qm = 2.1561373280078153
Qt = 0.12985582342677932
end
ρ_b = [1050, 700] # % kg/m^3
d_b = [d_b_1, d_b_2] * 1e-6 # m
d_s = [d_s_1, d_s_2, d_s_3] * 1e-6
@parameters begin
# FIXED VARIABLE FOR PSC
A = 7.065 * 10^2 # ; % m2
A_t = 555
μ = 0.000325 # ; % kg/m *s
g = 9.81 # ; % m/s
# FIXED VARIABLE FOR FT CELL
Aft = 163 # ; % m2
cd = 0.47 # drag coefficient
# ADJUSTABLE
ρ_s = 2650 # [2650, 2650, 2650] # % kg/m^3
ρ_w = 971.8 # % kg/m^3
d_b_1 = 21
d_b_2 = 240
d_s_1 = 7.5
d_s_2 = 160
d_s_3 = 380
# Tunable parameters
k_fines = 1.31
n_bf_1 = 6.05
n_bf_2 = 5.6
n_sf_1 = 7.3
n_bmu_1 = 5.3
n_bmu_2 = 4.2
n_smu_1 = 6.4
n_bt_1 = 6.5
n_st_1 = 7
n_st_2 = 6
n_st_3 = 5
λ_bf = 1.0906647723893508
λ_ff = 1.0909709890659598
λ_bm = 0.9505347443893531
λ_sm = 1.0616129443794053
λ_bt = 0.7502931250710971
λ_st = 1.0505261676728082
λ_fm = 0.8900741387772585
λ_ft = 1.060947106931769
end
# Define the equations
psv_eqs = [# WATER BALANCE
0 ~ 1 - (α_wf + (α_bf_1 + α_bf_2) + α_sf_1), # 1
0 ~ 1 - (α_wmu + (α_bmu_1 + α_bmu_2) + α_smu_1), # 2
0 ~ α_wm - (1 - (α_bm_1 + α_bm_2) - (α_sm_1 + α_sm_2 + α_sm_3)), # 3
0 ~ 1 - (α_wt + α_bt_1 + (α_st_1 + α_st_2 + α_st_3)), # 4
# MIXTURE DENSITY
0 ~ ρ_f - (ρ_w * α_wf + ρ_b[1] * α_bf_1 + ρ_s * α_sf_1), # 5 # ρ_b[2] * α_bf_2 + ρ_s[2] * α_sf_2 + ρ_s[3] * α_sf_3 = 0
0 ~ ρ_mu - (ρ_w * α_wmu + ρ_b[1] * α_bmu_1 + ρ_b[2] * α_bmu_2 + ρ_s * α_smu_1), # 6 # ρ_s * (α_smu_1 + α_smu_2 + α_smu_3)
0 ~ ρ_m - (ρ_w * α_wm + ρ_b[1] * α_bm_1 + ρ_b[2] * α_bm_2 + ρ_s * (α_sm_1 + α_sm_2 + α_sm_3)), # 7 ρ_s * (α_sm_1 + α_sm_2 + α_sm_3)
0 ~ ρ_t - (ρ_w * α_wt + ρ_b[1] * α_bt_1 + ρ_s * (α_st_1 + α_st_2 + α_st_3)), # 8 ρ_b[2] * α_bf_2
# Water (Bulf flow) Velocity
0 ~ -v_wf - (Qf + A * (α_bmu_1 * v_bf_1 + α_bmu_2 * v_bf_2 + α_smu_1 * v_sf_1)) / (α_wf * A) , # 9 + α_smu_2 * v_sf_2 + α_smu_3 * v_sf_3
0 ~ -v_wmu - (Qf + A * (α_bm_1 * v_bmu_1 + α_bm_2 * v_bmu_2 + α_sm_1 * v_smu_1)) / (α_wmu * A), # 10 + α_sm_2 * v_smu_2 + α_sm_3 * v_smu_3
0 ~ -v_wt + (Qt - A_t * (α_bm_1 * v_bt_1 + α_sm_1 * v_st_1 + α_sm_2 * v_st_2 + α_sm_3 * v_st_3)) / (α_wt * A_t), # 11 + α_bm_2 * v_bt_2
######### SETTLING VELOCITY ##########
# Froth
0 ~ v_bf_1 - (g * (d_b[1] ^ 2) * (ρ_b[1] - ρ_f) / (18 * μ) * max((α_wf - (k_fines - 1) * α_sf_1), 0) ^ n_bf_1 / α_wf ^ 2 + v_wf), # 12
0 ~ v_bf_2 - (g * (d_b[2] ^ 2) * (ρ_b[2] - ρ_f) / (18 * μ) * max((α_wf - (k_fines - 1) * α_sf_1), 0) ^ n_bf_2 / α_wf ^ 2 + v_wf), # 13
0 ~ v_sf_1 - (g * (d_s[1] ^ 2) * (ρ_s - ρ_f) / (18 * μ) * max((α_wf - (k_fines - 1) * α_sf_1), 0) ^ n_sf_1 / α_wf ^ 2 + v_wf), # 14
# MIDDLING UPPER
0 ~ v_bmu_1 - (g * (d_b[1] ^ 2) * (ρ_b[1] - ρ_mu) / (18 * μ) * max((α_wmu - (k_fines - 1) * α_smu_1), 0) ^ n_bmu_1 / α_wmu ^ 2 + v_wmu), # 15
0 ~ v_bmu_2 - (g * (d_b[2] ^ 2) * (ρ_b[2] - ρ_mu) / (18 * μ) * max((α_wmu - (k_fines - 1) * α_smu_1), 0) ^ n_bmu_2 / α_wmu ^ 2 + v_wmu), # 16
0 ~ v_smu_1 - (g * (d_s[1] ^ 2) * (ρ_s - ρ_mu) / (18 * μ) * max((α_wmu - (k_fines - 1) * α_smu_1), 0) ^ n_smu_1 / α_wmu ^ 2 + v_wmu), # 17
# Tailing
0 ~ v_bt_1 - (g * (d_b[1] ^ 2) * (ρ_b[1] - ρ_t) / (18 * μ) * max((α_wt - (k_fines - 1) * α_st_1), 0) ^ n_bt_1 / α_wt ^ 2 + v_wt), # 18
0 ~ v_st_1 - (g * (d_s[1] ^ 2) * (ρ_s - ρ_t) / (18 * μ) * max((α_wt - (k_fines - 1) * α_st_1), 0) ^ n_st_1 / α_wt ^ 2 + v_wt), # 19
0 ~ v_st_2 - (g * (d_s[2] ^ 2) * (ρ_s - ρ_t) / (18 * μ) * max((α_wt - (k_fines - 1) * α_st_1), 0) ^ n_st_2 / α_wt ^ 2 + v_wt), # 20
0 ~ v_st_3 - (g * (d_s[3] ^ 2) * (ρ_s - ρ_t) / (18 * μ) * max((α_wt - (k_fines - 1) * α_st_1), 0) ^ n_st_3 / α_wt ^ 2 + v_wt), # 21
########## VOLUME FRACTION CHANGE ##########
# % Froth
0 ~ α_bf_1 + (α_bmu_1 * A * v_bf_1) / (λ_bf * Qf), # 22
0 ~ α_bf_2 + (α_bmu_2 * A * v_bf_2) / (λ_bf * Qf), # 23
0 ~ α_sf_1 + (α_smu_1 * A * v_sf_1) / (λ_ff * Qf), # 24
# % MIDDLING UPPER
0 ~ α_bmu_1 - (α_bm_1 * A * v_bmu_1) / (A * v_bf_1), # 25
0 ~ α_bmu_2 - (α_bm_2 * A * v_bmu_2) / (A * v_bf_2), # 26
0 ~ α_smu_1 - (α_sm_1 * A * v_smu_1) / (A * v_sf_1), # 27
# % MIDDLING (SOURCE ZONE)
0 ~ (Qfd * α_bfd_1 + Qff * 0 - λ_bm * Qm * α_bm_1 - α_bm_1 * A_t * v_bt_1 + α_bm_1 * A * v_bmu_1), # 28 α_bff_1
0 ~ (Qfd * α_bfd_2 + Qff * 0 - λ_bm * Qm * α_bm_2 - α_bm_2 * A_t * v_bt_2 + α_bm_2 * A * v_bmu_2), # 29 α_bff_2
0 ~ (Qfd * α_sfd_1 + Qff * 0 - λ_bm * Qm * α_sm_1 - α_sm_1 * A_t * v_st_1 + α_sm_1 * A * v_smu_1), # 30 α_sff_1
0 ~ (Qfd * α_sfd_2 + Qff * 0 - λ_bm * Qm * α_sm_2 - α_sm_2 * A_t * v_st_2 + α_sm_2 * A * v_smu_2), # 31 α_sff_2
0 ~ (Qfd * α_sfd_3 + Qff * 0 - λ_bm * Qm * α_sm_3 - α_sm_3 * A_t * v_st_3 + α_sm_3 * A * v_smu_3), # 32 α_sff_3
# TAILING
0 ~ α_st_1 - (α_sm_1 * A_t * v_st_1) / (λ_ft * Qt), # 33
0 ~ α_st_2 - (α_sm_2 * A_t * v_st_2) / (λ_st * Qt), # 34
0 ~ α_st_3 - (α_sm_3 * A_t * v_st_3) / (λ_st * Qt), # 35
0 ~ α_bt_1 - (α_bm_1 * A_t * v_bt_1) / (λ_bt * Qt), # 36
########### % COMBINED VOLUME FRACTION ##########
0 ~ α_bf + (α_bmu_1 * A * v_bf_1 + α_bmu_2 * A * v_bf_2) / (λ_bf * Qf), # 37
0 ~ α_sf - α_sf_1, # 38
0 ~ α_bmu - (α_bmu_1 + α_bmu_2), # 39
0 ~ α_smu - α_smu_1, # 40
0 ~ α_bm - (α_bm_1 + α_bm_2), # 41
0 ~ α_sm - (α_sm_1 + α_sm_2 + α_sm_3), # 42
0 ~ α_bt - α_bt_1, # 43
0 ~ α_st - (α_st_1 + α_st_2 + α_st_3), # 44
]
psv_state = [v_bf_1, v_bf_2, v_sf_1, v_bmu_1, v_bmu_2, v_smu_1, v_bt_1, v_st_1, v_st_2, v_st_3, v_wf, v_wmu, v_wt,
ρ_f, ρ_mu, ρ_m, ρ_t, α_bf, α_sf, α_bmu, α_smu, α_bm, α_sm, α_bt, α_st, α_bf_1, α_bf_2, α_sf_1,
α_bmu_1, α_bmu_2, α_smu_1, α_bm_1, α_bm_2, α_sm_1, α_sm_2, α_sm_3, α_bt_1, α_st_1, α_st_2, α_st_3, α_wf, α_wmu, α_wm, α_wt]
psv_param = [A, A_t, μ, g, ρ_s, ρ_w, d_b_1, d_b_2, d_s_1, d_s_2, d_s_3, k_fines, n_bf_1, n_bf_2, n_sf_1, n_bmu_1, n_bmu_2, n_smu_1, n_bt_1, n_st_1, n_st_2, n_st_3, λ_bf, λ_ff, λ_bm, λ_sm, λ_bt, λ_st, λ_fm, λ_ft]
@named psv_sys = NonlinearSystem(psv_eqs, psv_state, psv_param)
simplified_psv_sys = structural_simplify(psv_sys,; verbose=true)
But the number of equations and variables are equal in practice and everyting is defined properly, I guess.
What can I do to resolve the issue?
modelingtoolkit
nonlinearsolve