MATLAB: Model Migration from Code to Simulink: The Execution of Matlab Function Blocks

MATLABmigrationode solversimulation modelsimulink

I am in the process of migrating a Matlab model to SImulink. I do this by using Matlab's function block, as can be seen below:
For this block to work, I first execute an "initialise" m-file, which sends the required values as busses to Simulink's inputs. It works quite well – but when executing the Simulink model (code given below), I receive the following error:
I can find no reason for this error, unless the solver does not execute the function from top to bottom. Any views on rectifying this? The variable given as "undefined" is quite clearly defined and is not in an _if _statement.
The code (I know it is quite hefty, but I need only general comments):
function output=calculate_autoclave_odes(t,X,data,Q_removal,flow_constants)
%This function generates the set of ordinary differential equations that %need to be solved in order to calculate the dynamic behaviour of the %autoclave in the base metal refinery (i.e. the flow rates, compositions, %and temperatures of the respective streams)
residue_solid_fraction=zeros(1,7); spent_concentration = zeros(1,7); m_Liquid_in=zeros(1,8); m_solids_in=zeros(1,8); n_phases_in=zeros(1,14); m_phases_in=zeros(1,14); m_Liquid_TK50=zeros(1,8); P_water_2nd_stage = 0; Cp_L_9=0; extent_2 = zeros(1,21); %Initialise the vector containing the extents of reaction
%Assign input data to local variables %Composition of first stage solid residue residue_solid_fraction(1) = data.residue.Cu/100; %Mass fraction Cu in first stage residue residue_solid_fraction(2) = data.residue.Ni/100; %Mass fraction Ni in first stage residue residue_solid_fraction(3) = data.residue.Fe/100; %Mass fraction Fe in first stage residue residue_solid_fraction(4) = data.residue.S/100; %Mass fraction S in first stage residue residue_solid_fraction(5) = data.residue.Rh/1000000; %Mass fraction Rh in first stage residue residue_solid_fraction(6) = data.residue.Ru/1000000; %Mass fraction Ru in first stage residue residue_solid_fraction(7) = data.residue.Ir/1000000; %Mass fraction Ir in first stage residue Cu_as_Cu9S5 = data.residue.Cu_as_Cu9S5; %Fraction of copper in first stage residue present as Cu9S5, remainder as CuS Ni_as_NiS = data.residue.Ni_as_NiS; %Fraction of nickel in first stage residue present as NiS, remainder as Ni3S4 Rh_as_Rh2S3 = data.residue.Rh_as_Rh2S3; %Fraction of rhodium in first stage residue present as Rh2S3, remainder as Rh Ru_as_RuS2 = data.residue.Ru_as_RuS2; %Fraction of ruthenium in first stage residue present as RuS2, remainder as Ru Ir_as_Ir2S3 = data.residue.Ir_as_Ir2S3; %Fraction of iridium in first stage residue present as Ir2S3, remainder as Ir Fe_as_FeOHSO4 = data.residue.Fe_as_FeOHSO4; %Fraction of iron in first stage residue present as leachable Fe(OH)SO4, remainder not leachable iron species CuS_leachable = data.residue.CuS_leachable;
%Composition of spent electrolyte
spent_concentration(1) = data.spent.Cu; %Concentration of copper in spent electrolyte, g/l
spent_concentration(2) = data.spent.Ni; %Concentration of nickel in spent electrolyte, g/l
spent_concentration(3) = data.spent.Fe/1000; %Concentration of iron in spent electrolyte, g/l
spent_concentration(4) = data.spent.Rh/1000; %Concentration of rhodium in spent electrolyte, g/l
spent_concentration(5) = data.spent.Ru/1000; %Concentration of ruthenium in spent electrolyte, g/l
spent_concentration(6) = data.spent.Ir/1000; %Concentration of iridium in spent electrolyte, g/l
spent_concentration(7) = data.spent.Sulphuric_acid; %Concentration of sulphuric acid in spent electrolyte, g/l
%Component properties (densities, heat capacities, molecular weights etc.)
density_acid = data.acid.density; %Density of acid, kg/l
density_spent = data.spent.density; %Density of spent electrolyte, kg/l
density_residue = data.residue.density; %Density of first stage residue, kg/l
acid_concentration = data.acid.concentration; %Mass percentage H2SO4 in acid
MW = data.constants.MW; %Data set containing the molecular weights of individual elements, g/mol
Cp_H2SO4.a = 1000*data.constants.Cp.H2SO4_a/(2*MW.H+MW.S+4*MW.O); %Heat capacity constant of sulphuric acid, Cp = a+bT, T in dC, Cp in kJ/kg.dC

Cp_H2SO4.b = 1000*data.constants.Cp.H2SO4_b/(2*MW.H+MW.S+4*MW.O); %Heat capacity constant of sulphuric acid, Cp = a+bT, T in dC, Cp in kJ/kg.dC
Cp_H2O = 1000*data.constants.Cp.H2O/(2*MW.H+MW.O); %Heat capacity of water kJ/kg.dC
Cp_steam.a = 1000*data.constants.Cp.Steam_a/(2*MW.H+MW.O); %Heat capacity of steam, Cp = a+bT+cT^2+dT^3, T in dC, Cp in kJ/kg.dC



Cp_steam.b = 1000*data.constants.Cp.Steam_b/(2*MW.H+MW.O); %Heat capacity of steam, Cp = a+bT+cT^2+dT^3, T in dC, Cp in kJ/kg.dC
Cp_steam.c = 1000*data.constants.Cp.Steam_c/(2*MW.H+MW.O); %Heat capacity of steam, Cp = a+bT+cT^2+dT^3, T in dC, Cp in kJ/kg.dC
Cp_steam.d = 1000*data.constants.Cp.Steam_d/(2*MW.H+MW.O); %Heat capacity of steam, Cp = a+bT+cT^2+dT^3, T in dC, Cp in kJ/kg.dC
Cp_O2.a = 1000*data.constants.Cp.Oxygen_a/(2*MW.O); %Heat capacity of oxygen, Cp = a+bT+cT^2+dT^3, T in dC, Cp in kJ/kg.dC
Cp_O2.b = 1000*data.constants.Cp.Oxygen_b/(2*MW.O); %Heat capacity of oxygen, Cp = a+bT+cT^2+dT^3, T in dC, Cp in kJ/kg.dC


Cp_O2.c = 1000*data.constants.Cp.Oxygen_c/(2*MW.O); %Heat capacity of oxygen, Cp = a+bT+cT^2+dT^3, T in dC, Cp in kJ/kg.dC
Cp_O2.d = 1000*data.constants.Cp.Oxygen_d/(2*MW.O); %Heat capacity of oxygen, Cp = a+bT+cT^2+dT^3, T in dC, Cp in kJ/kg.dC
Cp_NiS = data.constants.Cp.NiS/(MW.Ni+MW.S); %Heat capacity of NiS, kJ/kg.K
Cp_CuS = data.constants.Cp.CuS/(MW.Cu+MW.S); %Heat capacity of CuS, kJ/kg.K
Cp_Cu9S5 = data.constants.Cp.Cu9S5/(9*MW.Cu+5*MW.S); %Heat capacity of Cu9S5, kJ/kg.K
Cp_Ni3S4.a = data.constants.Cp.Ni3S4_a/(3*MW.Ni+4*MW.S); %Heat capacity of Ni3S4, Cp = a+bT+cT^2+dT^3+e/(T^2), T in K, Cp in kJ/kg.K




Cp_Ni3S4.b = data.constants.Cp.Ni3S4_b/(3*MW.Ni+4*MW.S); %Heat capacity of Ni3S4, Cp = a+bT+cT^2+dT^3+e/(T^2), T in K, Cp in kJ/kg.K
Cp_Ni3S4.c = data.constants.Cp.Ni3S4_c/(3*MW.Ni+4*MW.S); %Heat capacity of Ni3S4, Cp = a+bT+cT^2+dT^3+e/(T^2), T in K, Cp in kJ/kg.K
Cp_Ni3S4.d = data.constants.Cp.Ni3S4_d/(3*MW.Ni+4*MW.S); %Heat capacity of Ni3S4, Cp = a+bT+cT^2+dT^3+e/(T^2), T in K, Cp in kJ/kg.K
Cp_Ni3S4.e = data.constants.Cp.Ni3S4_e/(3*MW.Ni+4*MW.S); %Heat capacity of Ni3S4, Cp = a+bT+cT^2+dT^3+e/(T^2), T in K, Cp in kJ/kg.K
%Heats of state changes (heats of reactions, heat of evaporation)
H_Evap_H2O = data.constants.H_Evap_H2O; %Heat of vaporisation for water at 100dC, kJ/kg
H_Reaction = data.constants.H_reaction; %Data set containing the heat of reaction for the respective reactions, kJ/mol
%Stream information (flow rates, temperatures)
m_S_1 = data.flow.stream_1.mass_flow_solids; %Mass flow rate of solids in stream 1, kg/h
V_spent_1 = data.flow.stream_1.volume_flow_spent; %Volumetric flow rate of spent electrolyte in stream 1, l/h
V_2 = data.flow.stream_2.volume_flow; %Volumetric flow rate of stream 2 (spent electrolyte), l/h
V_3 = data.flow.stream_3.volume_flow; %Volumetric flow rate of stream 3 (water), l/h
V_4 = data.flow.stream_4.volume_flow; %Volumetric flow rate of stream 4 (acid), l/h
m_9 = data.flow.stream_9.mass_flow; %Mass flow rate of stream 9, kg/h
T_5 = data.flow.stream_5.temperature; %Temperature of the stream 5, dC
T_21 = data.flow.stream_21.temperature;
fraction_oxygen_excess = data.flow.oxygen.total_excess; %Fraction excess oxygen fed to autoclave
fraction_oxygen_stream_10 = data.flow.oxygen.fraction_10; %Fraction of total oxygen fed in stream 10
fraction_oxygen_stream_11 = data.flow.oxygen.fraction_11; %Fraction of total oxygen fed in stream 11
V_18 = data.flow.stream_18.volume; %Volumetric flow rate of spent electrolyte in stream 18, l/h
V_19 = data.flow.stream_19.volume; %Volumetric flow rate of water in stream 19, l/h
V_20 = data.flow.stream_20.volume; %Volumetric flow rate of acid in stream 20, l/h
Q_Removed_2 = Q_removal.compartment_2; %Rate at which heat is removed from autoclave compartment 2, kJ/h
Q_Removed_3 = Q_removal.compartment_3; %Rate at which heat is removed from autoclave compartment 3,kJ/h
m_13 = Q_removal.m_13; %Mass flow rate of steam to compartment 3, kg/h
%Autoclave conditions (pressure, volume, etc.)
autoclave_P = data.equipment.autoclave_pressure; %Operating pressure of the autoclave, bar
oxygen_T = data.flow.oxygen_temperature; %Oxygen temperature, dC
steam_T = data.flow.steam_temperature; %Steam temperature, dC
volume_1 = 1000*data.equipment.autoclave_volume.compartment_1; %Volume of first autoclave compartment in l (second stage leach)
volume_2 = 1000*data.equipment.autoclave_volume.compartment_2; %Volume of second autoclave compartment in l (second stage leach)
volume_3 = 1000*data.equipment.autoclave_volume.compartment_3; %Volume of third autoclave compartment in l (second stage leach)
volume_4 = 1000*data.equipment.autoclave_volume.compartment_4; %Volume of fourth autoclave compartment in l (third stage leach)
agitator_power = 3600*data.equipment.agitator_power; %Power of agitator drive, kJ/h
autoclave_outer_diameter = data.equipment.autoclave.outer_diameter; %Outside diameter of autoclave, m
autoclave_surface_temperature = data.equipment.autoclave.outer_temperature; %Temperature of autoclave outer surface, dC
autoclave_emissivity = data.equipment.autoclave.emissivity; %Emissivity of autoclave outer surface
length_1 = data.equipment.autoclave.compartment_1_length; %Length of the first compartment, m
length_2 = data.equipment.autoclave.compartment_2_length; %Length of the second compartment, m
length_3 = data.equipment.autoclave.compartment_3_length; %Length of the third compartment, m
length_4 = data.equipment.autoclave.compartment_4_length; %Length of the fourth compartment, m
air_temperature = data.constants.air.temperature; %Ambient temperature, dC
air_k = data.constants.air.k; %Thermal conductivity of air, W/m.K
air_viscosity = data.constants.air.viscosity; %Viscosity of air, m^2/s
air_thermal_diff = data.constants.air.thermal_diff; %Thermal diffusivity of air, m^2/s
air_Pr = data.constants.air.Pr; %Prandtl number of air at 300K
air_thermal_expansion = data.constants.air.thermal_expansion; %Thermal expansion coefficient, K^-1
%Estimate the constants to relate the outlet flow rate from a tank to the mass of slurry in the tank
prep_tank_constant = flow_constants.prep_tank_2;
recycle_tank_constant = flow_constants.recycle_tank;
discharge_tank_constant = flow_constants.discharge_tank;
thickener_constant = flow_constants.thickener;
prep_tank_3_constant = flow_constants.prep_tank_3;
%Calculation parameters
max_loops = 12000; %Maximum number of iterative calculation loops allowed
tolerance_percent = 0.0001; %Maximum percentage by which calculated temperature of composition may vary from initial guess
%Specify step changes in process conditions at specific time instances by varying any of the operating variables above (THIS MUST BE REPLICATED IN THE FUNCTION 'calculate_dynamic_behaviour')
% if t<6; % V_4 = V_4; % elseif t<72; % V_4 = 1.5*V_4; % else % V_4 = V_4; % end;
% if t==0; % m_13 = m_13; % else % m_13 = m_13+10; % end;
for i=(1:length(X)) if (X(i)<0) X(i) =0; end; end;
%Calculate flow rates of streams fed to second stage slurry preparation tank, TK-10 m_spent_1 = V_spent_1*density_spent; %Mass flow rate of spent electrolyte in stream 1, kg/h m_2 = V_2*density_spent; %Mass flow rate of spent electrolyte in stream 2, kg/h m_3 = V_3; %Mass flow rate of water in stream 3, kg/h m_4 = V_4*density_acid; %Mass flow rate of stream 4 (acid), kg/h m_H2SO4_4 = (acid_concentration/100)*m_4; %Mass flow rate of pure H2SO4 in stream 4, kg/h
%Calculate the total flow rates of liquids and dissolved species entering the second stage slurry preparation tank, TK-10, kg/h m_Liquid_in(1) = m_3 + m_2 + m_spent_1 + m_4; %Mass balance to calculate total flow rate of aqueous part entering TK-10, kg/h m_Liquid_in(2:7) = (V_2+V_spent_1)*spent_concentration(1:6)/1000; %Mass balance to calculate flow rate of dissolved metals entering TK-10, kg/h m_Liquid_in(8) = m_H2SO4_4 + (V_2+V_spent_1)*spent_concentration(7)/1000; %Mass balance to calculate flow rate of acid entering TK-10, kg/h
%Calculate the flow rates of solid elemental species entering the second stage slurry preparation tank, TK-10, kg/h m_solids_in(1) = m_S_1; %Total solids into preparation tank, TK-10, kg/h m_solids_in(2:8) = m_solids_in(1)*residue_solid_fraction; %Calculate flow rate of individual elements;composition of solids assumed to remain unchanged in the preparation tank 400TK-010, kg/h
%Calculate the molar amounts of the solid phases entering the second stage slurry preparation tank, TK-10 (NiS, Ni3S4, Cu9S5, CuS, Fe(OH)SO4, Rh2S3, Rh, RhO2, RuS2, Ru, RuO2, Ir2S3, Ir, IrO2) n_phases_in(1) = (m_solids_in(3)/MW.Ni)*Ni_as_NiS; %Molar flow rate of NiS in stream 5, kmol/h n_phases_in(2) = ((m_solids_in(3)/MW.Ni)*(1-Ni_as_NiS))/3; %Molar flow rate of Ni3S4 in stream 5, kmol/h n_phases_in(3) = ((m_solids_in(2)/MW.Cu)*Cu_as_Cu9S5)/9; %Molar flow rate of Cu9S5 in stream 5, kmol/h n_phases_in(4) = (m_solids_in(2)/MW.Cu)*(1-Cu_as_Cu9S5); %Molar flow rate of CuS in stream 5, kmol/h n_phases_in(5) = (m_solids_in(4)/MW.Fe)*Fe_as_FeOHSO4; %Molar flow rate of FeOHSO4 in stream 5, kmol/h n_phases_in(6) = ((m_solids_in(6)/MW.Rh)*Rh_as_Rh2S3)/2; %Molar flow rate of Rh2S3 in stream 5, kmol/h n_phases_in(7) = (m_solids_in(6)/MW.Rh)*(1-Rh_as_Rh2S3); %Molar flow rate of Rh in stream 5, kmol/h n_phases_in(8) = 0; %Molar flow rate of RhO2 in stream 5, kmol/h n_phases_in(9) = ((m_solids_in(7)/MW.Ru)*Ru_as_RuS2); %Molar flow rate of RuS2 in stream 5, kmol/h n_phases_in(10) = (m_solids_in(7)/MW.Ru)*(1-Ru_as_RuS2); %Molar flow rate of Ru in stream 5, kmol/h n_phases_in(11) = 0; %Molar flow rate of RuO2 in stream 5, kmol/h n_phases_in(12) = ((m_solids_in(8)/MW.Ir)*Ir_as_Ir2S3)/2; %Molar flow rate of Ir2S3 in stream 5, kmol/h n_phases_in(13) = (m_solids_in(8)/MW.Ir)*(1-Ir_as_Ir2S3); %Molar flow rate of Ir in stream 5, kmol/h n_phases_in(14) = 0; %Molar flow rate of IrO2 in stream 5, kmol/h
%Calculate the mass amounts of the solid phases entering the second stage slurry preparation tank, TK-10 (NiS, Ni3S4, Cu9S5, CuS, Fe(OH)SO4, Rh2S3, Rh, RhO2, RuS2, Ru, RuO2, Ir2S3, Ir, IrO2) m_phases_in(1) = n_phases_in(1)*(MW.Ni+MW.S); %Mass flow rate of NiS in stream 5, kg/h m_phases_in(2) = n_phases_in(2)*(3*MW.Ni+4*MW.S); %Mass flow rate of Ni3S4 in stream 5, kg/h m_phases_in(3) = n_phases_in(3)*(9*MW.Cu+5*MW.S); %Mass flow rate of Cu9S5 in stream 5, kg/h m_phases_in(4) = n_phases_in(4)*(MW.Cu+MW.S); %Mass flow rate of CuS in stream 5, kg/h m_phases_in(5) = n_phases_in(5)*(MW.Fe+5*MW.O+MW.H+MW.S); %Mass flow rate of FeOHSO4 in stream 5, kg/h m_phases_in(6) = n_phases_in(6)*(2*MW.Rh+3*MW.S); %Mass flow rate of Rh2S3 in stream 5, kg/h m_phases_in(7) = n_phases_in(7)*(MW.Rh); %Mass flow rate of Rh in stream 5, kg/h m_phases_in(8) = n_phases_in(8)*(MW.Rh+2*MW.O); %Mass flow rate of RhO2 in stream 5, kg/h m_phases_in(9) = n_phases_in(9)*(MW.Ru+2*MW.S); %Mass flow rate of RuS2 in stream 5, kg/h m_phases_in(10) = n_phases_in(10)*(MW.Ru); %Mass flow rate of Ru in stream 5, kg/h m_phases_in(11) = n_phases_in(11)*(MW.Ru+2*MW.O); %Mass flow rate of RuO2 in stream 5, kg/h m_phases_in(12) = n_phases_in(12)*(2*MW.Ir+3*MW.S); %Mass flow rate of Ir2S3 in stream 5, kg/h m_phases_in(13) = n_phases_in(13)*(MW.Ir); %Mass flow rate of Ir in stream 5, kg/h m_phases_in(14) = n_phases_in(14)*(MW.Ir+2*MW.O); %Mass flow rate of IrO2 in stream 5, kg/h
%Calculate the total flow rates of liquids and dissolved species entering the third stage slurry preparation tank, TK-50, kg/h m_18 = V_18*density_spent; %Mass flow rate of spent electrolyte in stream 2, kg/h m_19 = V_19; %Mass flow rate of water in stream 3, kg/h m_20 = V_20*density_acid; %Mass flow rate of stream 4 (acid), kg/h m_H2SO4_20 = (acid_concentration/100)*m_20; %Mass flow rate of pure H2SO4 in stream 4, kg/h
%Calculate the total flow rates of liquids and dissolved species entering the third stage slurry preparation tank, TK-50, kg/h m_Liquid_TK50(1) = m_18 + m_19 + m_20; %Mass balance to calculate total flow rate of aqueous part entering TK-50, kg/h m_Liquid_TK50(2:7) = (V_18)*spent_concentration(1:6)/1000; %Mass balance to calculate flow rate of dissolved metals entering TK-50, kg/h m_Liquid_TK50(8) = m_H2SO4_20 + (V_18)*spent_concentration(7)/1000; %Mass balance to calculate flow rate of acid entering TK-50, kg/h
%Calculate the fraction oxygen in the vapour space of the autoclave (assume that the vapour phases in the different compartments have the same temperature, equal to the temperature in the first compartment) %Calculate water vapour pressure (total pressure equal to water vapour pressure plus oxygen partial pressure, assuming that vapour phase is saturated with water vapour) if X(73)<60 %Temperatures less than 60dC require different constants for Antoine equation P_water_2nd_stage = ((10^(8.10765-1750.286/(X(73)+235)))/760)*1.01325; %Vapour pressure of water, bar else P_water_2nd_stage = ((10^(7.96681-1668.210/(X(73)+228)))/760)*1.01325; %Vapour pressure of water,bar end; %Calculate mol fraction oxygen in vapour y_O2_2nd_stage = 1-P_water_2nd_stage/autoclave_P; %Mole fraction of oxygen in vapour of second stage leach, assuming that vapour is saturated with water vapour
%Estimate heat losses from respective autoclave compartments due to convection and radiation Ra = (9.81*air_thermal_expansion*(autoclave_surface_temperature-air_temperature)*(autoclave_outer_diameter^3))/(air_viscosity*air_thermal_diff); Nu = (0.6+(0.387*Ra^(1/6))/((1+(0.559/air_Pr)^(9/16))^(8/27)))^2; h = air_k*Nu/autoclave_outer_diameter; Q_loss = (h*22/7*autoclave_outer_diameter*(autoclave_surface_temperature-air_temperature)+autoclave_emissivity*22/7*autoclave_outer_diameter*5.67*(10^-8)*((autoclave_surface_temperature+273.15)^4-(air_temperature+273.15)^4))*3600/1000;%Heat loss per meter of autoclave, kJ/h.m
%Calculate the theoretical amount of oxygen required for complete dissolution of NiS and Cu9S5 in feed to process (assuming all copper present as Cu9S5 and all nickel present at NiS) n_Cu_solids_in = m_solids_in(2)/MW.Cu; %Molar flow rate of solid Cu in stream 5, kmol/h n_Ni_solids_in = m_solids_in(3)/MW.Ni; %Molar flow rate of solid Ni in stream 5, kmol/h theoretical_oxygen = 2*n_Ni_solids_in + 12*(n_Cu_solids_in)/9; %Theoretical amount of oxygen required for complete leaching of NiS and Cu9S5, kmol/h
%Calculate actual flow rates of the oxygen feed streams n_10 = theoretical_oxygen*fraction_oxygen_excess*fraction_oxygen_stream_10; %Flow rate of stream 10, kmol/h n_11 = theoretical_oxygen*fraction_oxygen_excess*fraction_oxygen_stream_11; %Flow rate of stream 11, kmol/h n_12 = theoretical_oxygen*fraction_oxygen_excess-n_10-n_11; %Flow rate of stream 12, kmol/h m_10 = n_10*(2*MW.O); %Oxygen mass flow rate in stream 10, kg/h m_11 = n_11*(2*MW.O); %Oxygen mass flow rate in stream 11, kg/h m_12 = n_12*(2*MW.O); %Oxygen mass flow rate in stream 12, kg/h
%Estimate the heat capacities of the liquid phase in the respective streams, kJ/kg.dC Cp_L_5 = X(24)*Cp_H2SO4.a + (1-X(24))*Cp_H2O; %Stream 5 Cp_L_7 = X(48)*Cp_H2SO4.a+(1-X(48))*Cp_H2O; %Stream 7 Cp_L_9 = X(72)*Cp_H2SO4.a+(1-X(72))*Cp_H2O; %Stream 9 Cp_L_AC1 = Cp_L_9; %Stream AC1 Cp_L_AC2 = X(96)*Cp_H2SO4.a+(1-X(96))*Cp_H2O; %Stream AC2 Cp_L_21 = X(193)*Cp_H2SO4.a+(1-X(193))*Cp_H2O; %Stream 21 Cp_L_22 = X(216)*Cp_H2SO4.a+(1-X(216))*Cp_H2O; %Stream 22 Cp_L_14 = X(120)*Cp_H2SO4.a+(1-X(120))*Cp_H2O; %Stream 14
%Normalise mass fractions of NiS, Ni3S4, CuS, and Cu9S5 in solid phase in the respective streams Normalx_NiS_S_5 = X(3)/(X(3)+X(4)+X(5)+X(6)); Normalx_Ni3S4_S_5 = X(4)/(X(3)+X(4)+X(5)+X(6)); Normalx_Cu9S5_S_5 = X(5)/(X(3)+X(4)+X(5)+X(6)); Normalx_CuS_S_5 = X(6)/(X(3)+X(4)+X(5)+X(6));
Normalx_NiS_S_7 = X(27)/(X(27)+X(28)+X(29)+X(30));
Normalx_Ni3S4_S_7 = X(28)/(X(27)+X(28)+X(29)+X(30));
Normalx_Cu9S5_S_7 = X(29)/(X(27)+X(28)+X(29)+X(30));
Normalx_CuS_S_7 = X(30)/(X(27)+X(28)+X(29)+X(30));
Normalx_NiS_S_9 = X(51)/(X(51)+X(52)+X(53)+X(54));
Normalx_Ni3S4_S_9 = X(52)/(X(51)+X(52)+X(53)+X(54));
Normalx_Cu9S5_S_9 = X(53)/(X(51)+X(52)+X(53)+X(54));
Normalx_CuS_S_9 = X(54)/(X(51)+X(52)+X(53)+X(54));
Normalx_NiS_S_AC2 = X(75)/(X(75)+X(76)+X(77)+X(78));
Normalx_Ni3S4_S_AC2 = X(76)/(X(75)+X(76)+X(77)+X(78));
Normalx_Cu9S5_S_AC2 = X(77)/(X(75)+X(76)+X(77)+X(78));
Normalx_CuS_S_AC2 = X(78)/(X(75)+X(76)+X(77)+X(78));
Normalx_NiS_S_21 = X(172)/(X(172)+X(173)+X(174)+X(175));
Normalx_Ni3S4_S_21 = X(173)/(X(172)+X(173)+X(174)+X(175));
Normalx_Cu9S5_S_21 = X(174)/(X(172)+X(173)+X(174)+X(175));
Normalx_CuS_S_21 = X(175)/(X(172)+X(173)+X(174)+X(175));
Normalx_NiS_S_22 = X(195)/(X(195)+X(196)+X(197)+X(198));
Normalx_Ni3S4_S_22 = X(196)/(X(195)+X(196)+X(197)+X(198));
Normalx_Cu9S5_S_22 = X(197)/(X(195)+X(196)+X(197)+X(198));
Normalx_CuS_S_22 = X(198)/(X(195)+X(196)+X(197)+X(198));
Normalx_NiS_S_14 = X(99)/(X(99)+X(100)+X(101)+X(102));
Normalx_Ni3S4_S_14 = X(100)/(X(99)+X(100)+X(101)+X(102));
Normalx_Cu9S5_S_14 = X(101)/(X(99)+X(100)+X(101)+X(102));
Normalx_CuS_S_14 = X(102)/(X(99)+X(100)+X(101)+X(102));
%Estimate heat capacity of the solid portion of the respective streams, kJ/kg.K Cp_S_5 = Normalx_NiS_S_5*Cp_NiS + Normalx_Cu9S5_S_5*Cp_Cu9S5 + Normalx_CuS_S_5*Cp_CuS + Normalx_Ni3S4_S_5*Cp_Ni3S4.a; Cp_S_7 = Normalx_NiS_S_7*Cp_NiS + Normalx_Cu9S5_S_7*Cp_Cu9S5 + Normalx_CuS_S_7*Cp_CuS + Normalx_Ni3S4_S_7*Cp_Ni3S4.a; Cp_S_9 = Normalx_NiS_S_9*Cp_NiS + Normalx_Cu9S5_S_9*Cp_Cu9S5 + Normalx_CuS_S_9*Cp_CuS + Normalx_Ni3S4_S_9*Cp_Ni3S4.a; Cp_S_AC1 = Cp_S_9; Cp_S_AC2 = Normalx_NiS_S_AC2*Cp_NiS + Normalx_Cu9S5_S_AC2*Cp_Cu9S5 + Normalx_CuS_S_AC2*Cp_CuS + Normalx_Ni3S4_S_AC2*Cp_Ni3S4.a; Cp_S_21 = Normalx_NiS_S_21*Cp_NiS + Normalx_Cu9S5_S_21*Cp_Cu9S5 + Normalx_CuS_S_21*Cp_CuS + Normalx_Ni3S4_S_21*Cp_Ni3S4.a; Cp_S_22 = Normalx_NiS_S_22*Cp_NiS + Normalx_Cu9S5_S_22*Cp_Cu9S5 + Normalx_CuS_S_22*Cp_CuS + Normalx_Ni3S4_S_22*Cp_Ni3S4.a; Cp_S_14 = Normalx_NiS_S_14*Cp_NiS + Normalx_Cu9S5_S_14*Cp_Cu9S5 + Normalx_CuS_S_14*Cp_CuS + Normalx_Ni3S4_S_14*Cp_Ni3S4.a;
%Calculate energy removal to cool liquid portion of stream 9 to 100dC, kJ/h Q_removed_L_9 = (m_9*X(65))*Cp_L_9*(X(73)-100);
%Calculate energy removal to cool solid portion of stream 9 to 100dC, kJ/h Q_removed_S_9 = (m_9*X(50))*Cp_S_9*(X(73)-100);
%Calculate total energy removal to cool stream 9 to 100dC, kJ/h Q_removal_9 = Q_removed_L_9+Q_removed_S_9;
%Calculate rate of water evaporation from stream 9 upon entering the flash recycle tank, kg/h %If recycle stream has a temperature less than 100dC, no evaporation will occur if (Q_removal_9>0) m_H2O_flash_evaporated = Q_removal_9/H_Evap_H2O; else m_H2O_flash_evaporated = 0; end; %If the recycle flow rate is too low (or the autoclave operating temperating temperature too high), all the water in the recycle stream will evaporate if (m_H2O_flash_evaporated>m_9*X(65)*(1-X(66)*(MW.Cu+MW.S+4*MW.O)/MW.Cu-X(67)*(MW.Ni+MW.S+4*MW.O)/MW.Ni-X(72))) m_H2O_flash_evaporated = m_9*X(65)*(1-X(66)*(MW.Cu+MW.S+4*MW.O)/MW.Cu-X(67)*(MW.Ni+MW.S+4*MW.O)/MW.Ni-X(72)); temp_flag = 1; %Raise indication flag that all the water in the recycle stream evaporates during flashing under these conditions else temp_flag = 0; end;
%Compartment 1 reaction rates and extents of reactions

%Estimate volumetric flow rate of the feed to the autoclave, stream 5
V_S_5 = X(2)*prep_tank_constant*X(1)^0.5/density_residue; %Volumetric flow rate of solids in stream 5, l/h
V_L_5 = X(17)*prep_tank_constant*X(1)^0.5/density_spent; %Volumetric flow rate of liquid in stream 5, l/h
V_5 = V_S_5 + V_L_5; %Total volumetric flow rate of stream 5, l/h
%Estimate volumetric flow rate of the feed to the autoclave, stream 7
V_S_7 = X(26)*recycle_tank_constant*X(25)^0.5/density_residue; %Volumetric flow rate of solids in stream 7, l/h
V_L_7 = X(41)*recycle_tank_constant*X(25)^0.5/density_spent; %Volumetric flow rate of liquid in stream 7, l/h
V_7 = V_S_7 + V_L_7; %Total volumetric flow rate of stream 7, l/h
%Estimate the volumetric flow rate of stream 9
V_S_9 = m_9*X(50)/density_residue; %Volumetric flow rate of solids in stream 9, l/h
V_L_9 = m_9*X(65)/density_spent; %Volumetric flow rate of liquid in stream 9, l/h
V_9 = V_S_9 + V_L_9; %Total volumetric flow rate of stream 9, l/h
%Calculate the concentrations of reagents present in compartment 1 in order to calculate the reaction rates
%Dissolved/liquid species
liq_fractions_1.H2SO4 = X(72); %Concentration of aqueous species in autoclave compartment 1, kg H2SO4/kg liquid
liq_fractions_1.Rh = X(69); %Concentration of aqueous species in autoclave compartment 1, kg Rh/kg liquid
liq_fractions_1.Ru = X(70); %Concentration of aqueous species in autoclave compartment 1, kg Ru/kg liquid
liq_fractions_1.Ir = X(71); %Concentration of aqueous species in autoclave compartment 1, kg Ir/kg liquid
%Solid species entering compartment 1
solids_in_1.NiS = X(26)*recycle_tank_constant*X(25)^0.5*X(27)/V_7; %Concentration of solids entering the first autoclave compartment, kg NiS/l slurry
solids_in_1.Ni3S4 = X(26)*recycle_tank_constant*X(25)^0.5*X(28)/V_7; %Concentration of solids entering the first autoclave compartment, kg Ni3S4/l slurry
solids_in_1.Cu9S5 = X(26)*recycle_tank_constant*X(25)^0.5*X(29)/V_7; %Concentration of solids entering the first autoclave compartment, kg Cu9S5/l slurry
solids_in_1.FeOHSO4 = X(26)*recycle_tank_constant*X(25)^0.5*X(31)/V_7; %Concentration of solids entering the first autoclave compartment, kg FeOHSO4/l slurry
solids_in_1.Rh2S3 = X(26)*recycle_tank_constant*X(25)^0.5*X(32)/V_7; %Concentration of solids entering the first autoclave compartment, kg Rh2S3/l slurry
solids_in_1.Rh = X(26)*recycle_tank_constant*X(25)^0.5*X(33)/V_7; %Concentration of solids entering the first autoclave compartment, kg Rh/l slurry
solids_in_1.RuS2 = X(26)*recycle_tank_constant*X(25)^0.5*X(35)/V_7; %Concentration of solids entering the first autoclave compartment, kg RuS2/l slurry
solids_in_1.Ru = X(26)*recycle_tank_constant*X(25)^0.5*X(36)/V_7; %Concentration of solids entering the first autoclave compartment, kg Ru/l slurry
solids_in_1.Ir2S3 = X(26)*recycle_tank_constant*X(25)^0.5*X(38)/V_7; %Concentration of solids entering the first autoclave compartment, kg Ir2S3/l slurry
solids_in_1.Ir = X(26)*recycle_tank_constant*X(25)^0.5*X(39)/V_7; %Concentration of solids entering the first autoclave compartment, kg Ir/l slurry
Cu9S5_fresh = X(2)*prep_tank_constant*X(1)^0.5*X(5)/V_5; %Concentration of Cu9S5 in the fresh feed to the flash recycle tank, kg Cu9S5/l slurry
%Solid species leaving compartment 1
solids_out_1.NiS = m_9*X(50)*X(51)/V_9; %Concentration of solids leaving the first autoclave compartment, kg NiS/l slurry
solids_out_1.Ni3S4 = m_9*X(50)*X(52)/V_9; %Concentration of solids leaving the first autoclave compartment, kg Ni3S4/l slurry
solids_out_1.Cu9S5 = m_9*X(50)*X(53)/V_9; %Concentration of solids leaving the first autoclave compartment, kg Cu9S5/l slurry
solids_out_1.CuS = m_9*X(50)*X(54)/V_9; %Concentration of solids leaving the first autoclave compartment, kg CuS/l slurry
solids_out_1.FeOHSO4 = m_9*X(50)*X(55)/V_9; %Concentration of solids leaving the first autoclave compartment, kg FeOHSO4/l slurry
solids_out_1.Rh2S3 = m_9*X(50)*X(56)/V_9; %Concentration of solids leaving the first autoclave compartment, kg Rh2S3/l slurry
solids_out_1.Rh = m_9*X(50)*X(57)/V_9; %Concentration of solids leaving the first autoclave compartment, kg Rh/l slurry
solids_out_1.RuS2 = m_9*X(50)*X(59)/V_9; %Concentration of solids leaving the first autoclave compartment, kg RuS2/l slurry
solids_out_1.Ru = m_9*X(50)*X(60)/V_9; %Concentration of solids leaving the first autoclave compartment, kg Ru/l slurry
solids_out_1.Ir2S3 = m_9*X(50)*X(62)/V_9; %Concentration of solids leaving the first autoclave compartment, kg Ir2S3/l slurry
solids_out_1.Ir = m_9*X(50)*X(63)/V_9; %Concentration of solids leaving the first autoclave compartment, kg Ir/l slurry
solids_out_1.RhO2 = m_9*X(50)*X(58)/V_9; %Concentration of solids leaving the first autoclave compartment, kg RhO2/l slurry
solids_out_1.RuO2 = m_9*X(50)*X(61)/V_9; %Concentration of solids leaving the first autoclave compartment, kg RuO2/l slurry
solids_out_1.IrO2 = m_9*X(50)*X(64)/V_9; %Concentration of solids leaving the first autoclave compartment, kg IrO2/l slurry
%Oxygen solubility (mol/kg) in stream 9
Cu_molal_9 = 1000*(X(66)/(1-X(66)*(MW.Cu+MW.S+4*MW.O)/MW.Cu-X(67)*(MW.Ni+MW.S+4*MW.O)/MW.Ni-X(72)))/MW.Cu; %Molal concentration of copper in stream 9, mol Cu per kg solvent
Ni_molal_9 = 1000*(X(67)/(1-X(66)*(MW.Cu+MW.S+4*MW.O)/MW.Cu-X(67)*(MW.Ni+MW.S+4*MW.O)/MW.Ni-X(72)))/MW.Ni; %Molal concentration of nickel in stream 9, mol Ni per kg solvent
H2SO4_molal_9 = 1000*(X(72)/(1-X(66)*(MW.Cu+MW.S+4*MW.O)/MW.Cu-X(67)*(MW.Ni+MW.S+4*MW.O)/MW.Ni-X(72)))/(2*MW.H+MW.S+4*MW.O); %Molal concentration of H2SO4 in stream 9, mol H2SO4 per kg solvent
O2_solubility_9 = (2*MW.O/1000)*(calculate_oxygen_solubility(X(73),(100000*autoclave_P),Cu_molal_9,Ni_molal_9,H2SO4_molal_9)); %Oxygen solubility in stream 9, kg O2 per kg water
molal_O2_solubility_9 = (O2_solubility_9*1000)/(2*MW.O); %Oxygen solubility in stream 9, mol O2 per kg water
%Calculate the reaction rates in compartment 1, mol/l.min
r_1 = calculate_reaction_rates(X(73),l

Best Answer

Niel: Don't mean to discourage you, but please post your comments on this thread instead of sending personal messages in reply (we highly discourage contacting contributors personally on this forum since our participation is voluntary and not part of our job responsibilities. Thanks or understanding!)
---Message from Niel follows---
Thank you for pointing out that the errors I received were "compile time" errors. By pressing Cntl+D, I indeed received the same errors. The reason I thought that the problem lies with Simulink is that the same code ran perfectly in Matlab's normal code format. I did make a few changes to the code in order to incorporate it into a Matlab function block. Why would the same piece of code suddenly given errors in Simulink?
--------------------------------
When your code runs in the MATLAB Function block, it actually generates C code from the MATLAB code and generates a MEX-file from it. It is this MEX-file that is used for execution of the block. The compilation process for generating the C code is a lot stricter than how the MATLAB interpeter works, which is why you're seeing the failure only when using the code in the block. You should however see the same compilation error if you attempt to generate code from the function using MATLAB Coder (outside of Simulink).
Now, coming back to the original question - it looks like your code is too long to read legibly even if you paste it here completely. Are you absolutely sure that the variable that the error complains about is assigned outside of the "if" section as well? To be safe, you could try assigning it to a value at the top of your function, and then let it get reassigned later.
Related Question