-
Notifications
You must be signed in to change notification settings - Fork 6
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
Merge pull request #194 from EcoExtreML/fix_issue_99
Fix issue 99
- Loading branch information
Showing
20 changed files
with
323 additions
and
267 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Binary file not shown.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,82 @@ | ||
| function [RHS, AirMatrices, SAVE] = assembleCoefficientMatrices(AirMatrices, SoilVariables, Delt_t, P_g) | ||
| %{ | ||
| Assemble the coefficient matrices of Equation 4.32 STEMMUS Technical | ||
| Notes, page 44, for dry air equation. | ||
| %} | ||
|
|
||
| C1 = AirMatrices.C1; | ||
| C2 = AirMatrices.C2; | ||
| C3 = AirMatrices.C3; | ||
| C4 = AirMatrices.C4; | ||
| C4_a = AirMatrices.C4_a; | ||
| C5 = AirMatrices.C5; | ||
| C5_a = AirMatrices.C5_a; | ||
| C6 = AirMatrices.C6; | ||
| C7 = AirMatrices.C7; | ||
|
|
||
| ModelSettings = io.getModelSettings(); | ||
| n = ModelSettings.NN; | ||
|
|
||
| % Alias of SoilVariables | ||
| SV = SoilVariables; | ||
|
|
||
| if ModelSettings.Thmrlefc | ||
| RHS(1) = -C7(1) + (C3(1, 1) * P_g(1) + C3(1, 2) * P_g(2)) / Delt_t ... | ||
| - (C2(1, 1) / Delt_t + C5(1, 1)) * SV.TT(1) - (C2(1, 2) / Delt_t + C5(1, 2)) * SV.TT(2) ... | ||
| - (C1(1, 1) / Delt_t + C4(1, 1)) * SV.hh(1) - (C1(1, 2) / Delt_t + C4(1, 2)) * SV.hh(2) ... | ||
| + (C2(1, 1) / Delt_t) * SV.T(1) + (C2(1, 2) / Delt_t) * SV.T(2) ... | ||
| + (C1(1, 1) / Delt_t) * SV.h(1) + (C1(1, 2) / Delt_t) * SV.h(2); | ||
|
|
||
| for i = 2:ModelSettings.NL | ||
| ARG1 = C2(i - 1, 2) / Delt_t; | ||
| ARG2 = C2(i, 1) / Delt_t; | ||
| ARG3 = C2(i, 2) / Delt_t; | ||
|
|
||
| ARG4 = C1(i - 1, 2) / Delt_t; | ||
| ARG5 = C1(i, 1) / Delt_t; | ||
| ARG6 = C1(i, 2) / Delt_t; | ||
|
|
||
| RHS(i) = -C7(i) + (C3(i - 1, 2) * P_g(i - 1) + C3(i, 1) * P_g(i) + C3(i, 2) * P_g(i + 1)) / Delt_t ... | ||
| - (ARG1 + C5_a(i - 1)) * SV.TT(i - 1) - (ARG2 + C5(i, 1)) * SV.TT(i) - (ARG3 + C5(i, 2)) * SV.TT(i + 1) ... | ||
| - (ARG4 + C4_a(i - 1)) * SV.hh(i - 1) - (ARG5 + C4(i, 1)) * SV.hh(i) - (ARG6 + C4(i, 2)) * SV.hh(i + 1) ... | ||
| + ARG1 * SV.T(i - 1) + ARG2 * SV.T(i) + ARG3 * SV.T(i + 1) ... | ||
| + ARG4 * SV.h(i - 1) + ARG5 * SV.h(i) + ARG6 * SV.h(i + 1); | ||
| end | ||
|
|
||
| RHS(n) = -C7(n) + (C3(n - 1, 2) * P_g(n - 1) + C3(n, 1) * P_g(n)) / Delt_t ... | ||
| - (C2(n - 1, 2) / Delt_t + C5_a(n - 1)) * SV.TT(n - 1) - (C2(n, 1) / Delt_t + C5(n, 1)) * SV.TT(n) ... | ||
| - (C1(n - 1, 2) / Delt_t + C4_a(n - 1)) * SV.hh(n - 1) - (C1(n, 1) / Delt_t + C4(n, 1)) * SV.hh(n) ... | ||
| + (C2(n - 1, 2) / Delt_t) * SV.T(n - 1) + (C2(n, 1) / Delt_t) * SV.T(n) ... | ||
| + (C1(n - 1, 2) / Delt_t) * SV.h(n - 1) + (C1(n, 1) / Delt_t) * SV.h(n); | ||
| else | ||
| RHS(1) = -C7(1) + (C3(1, 1) * P_g(1) + C3(1, 2) * P_g(2)) / Delt_t ... | ||
| - (C1(1, 1) / Delt_t + C4(1, 1)) * SV.hh(1) - (C1(1, 2) / Delt_t + C4(1, 2)) * SV.hh(2) ... | ||
| + (C1(1, 1) / Delt_t) * SV.h(1) + (C1(1, 2) / Delt_t) * SV.h(2); | ||
| for i = 2:ModelSettings.NL | ||
| ARG4 = C1(i - 1, 2) / Delt_t; | ||
| ARG5 = C1(i, 1) / Delt_t; | ||
| ARG6 = C1(i, 2) / Delt_t; | ||
| RHS(i) = -C7(i) + (C3(i - 1, 2) * P_g(i - 1) + C3(i, 1) * P_g(i) + C3(i, 2) * P_g(i + 1)) / Delt_t ... | ||
| - (ARG4 + C4(i - 1, 2)) * SV.hh(i - 1) - (ARG5 + C4(i, 1)) * SV.hh(i) - (ARG6 + C4(i, 2)) * SV.hh(i + 1) ... | ||
| + ARG4 * SV.h(i - 1) + ARG5 * SV.h(i) + ARG6 * SV.h(i + 1); | ||
| end | ||
| RHS(n) = -C7(n) + (C3(n - 1, 2) * P_g(n - 1) + C3(n, 1) * P_g(n)) / Delt_t ... | ||
| - (C1(n - 1, 2) / Delt_t + C4(n - 1, 2)) * SV.hh(n - 1) - (C1(n, 1) / Delt_t + C4(n, 1)) * SV.hh(n) ... | ||
| + (C1(n - 1, 2) / Delt_t) * SV.h(n - 1) + (C1(n, 1) / Delt_t) * SV.h(n); | ||
| end | ||
|
|
||
| for i = 1:ModelSettings.NN | ||
| for j = 1:ModelSettings.nD | ||
| C6(i, j) = C3(i, j) / Delt_t + C6(i, j); | ||
| end | ||
| end | ||
|
|
||
| AirMatrices.C6 = C6; | ||
|
|
||
| SAVE(1, 1, 3) = RHS(1); | ||
| SAVE(1, 2, 3) = C6(1, 1); | ||
| SAVE(1, 3, 3) = C6(1, 2); | ||
| SAVE(2, 1, 3) = RHS(n); | ||
| SAVE(2, 2, 3) = C6(n - 1, 2); | ||
| SAVE(2, 3, 3) = C6(n, 1); | ||
| end |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,37 @@ | ||
| function [RHS, AirMatrices] = calculateBoundaryConditions(BoundaryCondition, AirMatrices, ForcingData, RHS, KT) | ||
| %{ | ||
| Determine the boundary condition for solving the dry air equation. | ||
| %} | ||
|
|
||
| TopPg = 100 .* (ForcingData.Pg_msr); | ||
| ModelSettings = io.getModelSettings(); | ||
| n = ModelSettings.NN; | ||
|
|
||
| % Apply the bottom boundary condition called for by NBCPB | ||
| if BoundaryCondition.NBCPB == 1 % Bounded bottom with the water table | ||
| RHS(1) = BtmPg; | ||
| AirMatrices.C6(1, 1) = 1; | ||
| RHS(2) = RHS(2) - AirMatrices.C6(1, 2) * RHS(1); | ||
| AirMatrices.C6(1, 2) = 0; | ||
| AirMatrices.C6_a(1) = 0; | ||
| elseif BoundaryCondition.NBCPB == 2 % The soil air is allowed to escape from the bottom | ||
| RHS(1) = RHS(1) + BoundaryCondition.BCPB; | ||
| end | ||
|
|
||
| % Apply the surface boundary condition called by NBCP | ||
| if BoundaryCondition.NBCP == 1 % Ponded infiltration with Bonded bottom | ||
| RHS(n) = BtmPg; | ||
| AirMatrices.C6(n, 1) = 1; | ||
| RHS(n - 1) = RHS(n - 1) - AirMatrices.C6(n - 1, 2) * RHS(n); | ||
| AirMatrices.C6(n - 1, 2) = 0; | ||
| AirMatrices.C6_a(n - 1) = 0; | ||
| elseif BoundaryCondition.NBCP == 2 % Specified flux on the surface | ||
| RHS(n) = RHS(n) - BoundaryCondition.BCP; | ||
| else | ||
| RHS(n) = TopPg(KT); | ||
| AirMatrices.C6(n, 1) = 1; | ||
| RHS(n - 1) = RHS(n - 1) - AirMatrices.C6(n - 1, 2) * RHS(n); | ||
| AirMatrices.C6(n - 1, 2) = 0; | ||
| AirMatrices.C6_a(n - 1) = 0; | ||
| end | ||
| end |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,73 @@ | ||
| function AirVariabes = calculateDryAirParameters(SoilVariables, GasDispersivity, TransportCoefficient, InitialValues, VaporVariables, ... | ||
| P_gg, Xah, XaT, Xaa, RHODA) | ||
| %{ | ||
| Calculate all the parameters related to dry air equation e.g., Equation | ||
| 3.59-3.64, STEMMUS Technical Notes, page 27-28. | ||
| %} | ||
| ModelSettings = io.getModelSettings(); | ||
| Constants = io.define_constants(); | ||
|
|
||
| AirVariabes.Cah = InitialValues.Cah; | ||
| AirVariabes.CaT = InitialValues.CaT; | ||
| AirVariabes.Caa = InitialValues.Caa; | ||
| AirVariabes.Kah = InitialValues.Kah; | ||
| AirVariabes.KaT = InitialValues.KaT; | ||
| AirVariabes.Kaa = InitialValues.Kaa; | ||
| AirVariabes.Vah = InitialValues.Vah; | ||
| AirVariabes.VaT = InitialValues.VaT; | ||
| AirVariabes.Vaa = InitialValues.Vaa; | ||
| AirVariabes.Cag = InitialValues.Cag; | ||
| AirVariabes.QL = InitialValues.QL; | ||
| AirVariabes.KLhBAR = InitialValues.KLhBAR; | ||
| AirVariabes.KLTBAR = InitialValues.KLTBAR; | ||
| AirVariabes.DhDZ = InitialValues.DhDZ; | ||
| AirVariabes.DTDZ = InitialValues.DTDZ; | ||
| GasDispersivity.DPgDZ = InitialValues.DPgDZ; | ||
|
|
||
| for i = 1:ModelSettings.NL | ||
| KLhBAR = (SoilVariables.KfL_h(i, 1) + SoilVariables.KfL_h(i, 2)) / 2; | ||
| KLTBAR = (InitialValues.KL_T(i, 1) + InitialValues.KL_T(i, 2)) / 2; | ||
| DDhDZ = (SoilVariables.hh(i + 1) - SoilVariables.hh(i)) / ModelSettings.DeltZ(i); | ||
| DhDZ = (SoilVariables.hh(i + 1) + SoilVariables.hh_frez(i + 1) - SoilVariables.hh(i) - SoilVariables.hh_frez(i)) / ModelSettings.DeltZ(i); | ||
| DTDZ = (SoilVariables.TT(i + 1) - SoilVariables.TT(i)) / ModelSettings.DeltZ(i); | ||
| DPgDZ = (P_gg(i + 1) - P_gg(i)) / ModelSettings.DeltZ(i); | ||
| DTDBAR = (TransportCoefficient.D_Ta(i, 1) + TransportCoefficient.D_Ta(i, 2)) / 2; | ||
|
|
||
| if SoilVariables.KLa_Switch == 1 | ||
| QL = -(KLhBAR * (DhDZ + DPgDZ / Constants.Gamma_w) + (KLTBAR + DTDBAR) * DTDZ + KLhBAR); | ||
| QL_h(i) = -(KLhBAR * (DhDZ + DPgDZ / Constants.Gamma_w) + KLhBAR); | ||
| QL_a(i) = -(KLhBAR * (DPgDZ / Constants.Gamma_w)); | ||
| QL_T(i) = -((KLTBAR + DTDBAR) * DTDZ); | ||
| else | ||
| QL = -(KLhBAR * DhDZ + (KLTBAR + DTDBAR) * DTDZ + KLhBAR); | ||
| QL_h(i) = -(KLhBAR * DhDZ + KLhBAR); | ||
| QL_T(i) = -((KLTBAR + DTDBAR) * DTDZ); | ||
| end | ||
|
|
||
| % used in EnrgyPARM | ||
| AirVariabes.KLhBAR(i) = KLhBAR; | ||
| AirVariabes.KLTBAR(i) = KLTBAR; | ||
| AirVariabes.DDhDZ(i) = DDhDZ; | ||
| AirVariabes.DhDZ(i) = DhDZ; | ||
| AirVariabes.DTDZ(i) = DTDZ; | ||
| AirVariabes.QL(i) = QL; | ||
|
|
||
| for j = 1:ModelSettings.nD | ||
| MN = i + j - 1; | ||
|
|
||
| AirVariabes.Cah(i, j) = Xah(MN) * (SoilVariables.POR(i) + (Constants.Hc - 1) * SoilVariables.Theta_LL(i, j)) + (Constants.Hc - 1) * RHODA(MN) * SoilVariables.DTheta_LLh(i, j); | ||
| AirVariabes.CaT(i, j) = XaT(MN) * (SoilVariables.POR(i) + (Constants.Hc - 1) * SoilVariables.Theta_LL(i, j)) + (Constants.Hc - 1) * RHODA(MN) * SoilVariables.DTheta_LLT(i, j); | ||
| AirVariabes.Caa(i, j) = Xaa(MN) * (SoilVariables.POR(i) + (Constants.Hc - 1) * SoilVariables.Theta_LL(i, j)); | ||
|
|
||
| AirVariabes.Kah(i, j) = Xah(MN) * (VaporVariables.D_V(i, j) + GasDispersivity.D_Vg(i)) + Constants.Hc * RHODA(MN) * SoilVariables.KfL_h(i, j); | ||
| AirVariabes.KaT(i, j) = XaT(MN) * (VaporVariables.D_V(i, j) + GasDispersivity.D_Vg(i)) + Constants.Hc * RHODA(MN) * (InitialValues.KL_T(i, j) + TransportCoefficient.D_Ta(i, j)); | ||
| AirVariabes.Kaa(i, j) = Xaa(MN) * (VaporVariables.D_V(i, j) + GasDispersivity.D_Vg(i)) + RHODA(MN) * (GasDispersivity.Beta_g(i, j) + Constants.Hc * SoilVariables.KfL_h(i, j) / Constants.Gamma_w); | ||
|
|
||
| AirVariabes.Cag(i, j) = Constants.Hc * RHODA(MN) * SoilVariables.KfL_h(i, j); | ||
|
|
||
| AirVariabes.Vah(i, j) = -(GasDispersivity.V_A(i) + Constants.Hc * QL / Constants.RHOL) * Xah(MN); | ||
| AirVariabes.VaT(i, j) = -(GasDispersivity.V_A(i) + Constants.Hc * QL / Constants.RHOL) * XaT(MN); | ||
| AirVariabes.Vaa(i, j) = -(GasDispersivity.V_A(i) + Constants.Hc * QL / Constants.RHOL) * Xaa(MN); | ||
| end | ||
| end | ||
| end |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,57 @@ | ||
| function AirMatrices = calculateMatricCoefficients(AirVariabes, InitialValues) | ||
| %{ | ||
| Calculate all the parameters related to matric coefficients e.g., | ||
| c1-c7 as in Equation 4.32 STEMMUS Technical Notes, page 44, which is | ||
| an example for soil moisture equation, but for dry air equation. | ||
| %} | ||
|
|
||
| ModelSettings = io.getModelSettings(); | ||
|
|
||
| AirMatrices.C1 = InitialValues.C1; | ||
| AirMatrices.C2 = InitialValues.C2; | ||
| AirMatrices.C3 = InitialValues.C3; | ||
| AirMatrices.C4 = InitialValues.C4; | ||
| AirMatrices.C5 = InitialValues.C5; | ||
| AirMatrices.C6 = InitialValues.C6; | ||
| AirMatrices.C7 = zeros(ModelSettings.NN); | ||
|
|
||
| for i = 1:ModelSettings.NL | ||
| AirMatrices.C1(i, 1) = AirMatrices.C1(i, 1) + AirVariabes.Cah(i, 1) * ModelSettings.DeltZ(i) / 2; | ||
| AirMatrices.C1(i + 1, 1) = AirMatrices.C1(i + 1, 1) + AirVariabes.Cah(i, 2) * ModelSettings.DeltZ(i) / 2; | ||
|
|
||
| AirMatrices.C2(i, 1) = AirMatrices.C2(i, 1) + AirVariabes.CaT(i, 1) * ModelSettings.DeltZ(i) / 2; | ||
| AirMatrices.C2(i + 1, 1) = AirMatrices.C2(i + 1, 1) + AirVariabes.CaT(i, 2) * ModelSettings.DeltZ(i) / 2; | ||
|
|
||
| AirMatrices.C3(i, 1) = AirMatrices.C3(i, 1) + AirVariabes.Caa(i, 1) * ModelSettings.DeltZ(i) / 2; | ||
| AirMatrices.C3(i + 1, 1) = AirMatrices.C3(i + 1, 1) + AirVariabes.Caa(i, 2) * ModelSettings.DeltZ(i) / 2; | ||
|
|
||
| C4ARG1 = (AirVariabes.Kah(i, 1) + AirVariabes.Kah(i, 2)) / (2 * ModelSettings.DeltZ(i)); | ||
| C4ARG2_1 = AirVariabes.Vah(i, 1) / 3 + AirVariabes.Vah(i, 2) / 6; | ||
| C4ARG2_2 = AirVariabes.Vah(i, 1) / 6 + AirVariabes.Vah(i, 2) / 3; | ||
|
|
||
| AirMatrices.C4(i, 1) = AirMatrices.C4(i, 1) + C4ARG1 - C4ARG2_1; | ||
| AirMatrices.C4(i, 2) = AirMatrices.C4(i, 2) - C4ARG1 - C4ARG2_2; | ||
| AirMatrices.C4(i + 1, 1) = AirMatrices.C4(i + 1, 1) + C4ARG1 + C4ARG2_2; | ||
| AirMatrices.C4_a(i) = -C4ARG1 + C4ARG2_1; | ||
|
|
||
| C5ARG1 = (AirVariabes.KaT(i, 1) + AirVariabes.KaT(i, 2)) / (2 * ModelSettings.DeltZ(i)); | ||
| C5ARG2_1 = AirVariabes.VaT(i, 1) / 3 + AirVariabes.VaT(i, 2) / 6; | ||
| C5ARG2_2 = AirVariabes.VaT(i, 1) / 6 + AirVariabes.VaT(i, 2) / 3; | ||
| AirMatrices.C5(i, 1) = AirMatrices.C5(i, 1) + C5ARG1 - C5ARG2_1; | ||
| AirMatrices.C5(i, 2) = AirMatrices.C5(i, 2) - C5ARG1 - C5ARG2_2; | ||
| AirMatrices.C5(i + 1, 1) = AirMatrices.C5(i + 1, 1) + C5ARG1 + C5ARG2_2; | ||
| AirMatrices.C5_a(i) = -C5ARG1 + C5ARG2_1; | ||
|
|
||
| C6ARG1 = (AirVariabes.Kaa(i, 1) + AirVariabes.Kaa(i, 2)) / (2 * ModelSettings.DeltZ(i)); | ||
| C6ARG2_1 = AirVariabes.Vaa(i, 1) / 3 + AirVariabes.Vaa(i, 2) / 6; | ||
| C6ARG2_2 = AirVariabes.Vaa(i, 1) / 6 + AirVariabes.Vaa(i, 2) / 3; | ||
| AirMatrices.C6(i, 1) = AirMatrices.C6(i, 1) + C6ARG1 - C6ARG2_1; | ||
| AirMatrices.C6(i, 2) = AirMatrices.C6(i, 2) - C6ARG1 - C6ARG2_2; | ||
| AirMatrices.C6(i + 1, 1) = AirMatrices.C6(i + 1, 1) + C6ARG1 + C6ARG2_2; | ||
| AirMatrices.C6_a(i) = -C6ARG1 + C6ARG2_1; | ||
|
|
||
| C7ARG = (AirVariabes.Cag(i, 1) + AirVariabes.Cag(i, 2)) / 2; | ||
| AirMatrices.C7(i) = AirMatrices.C7(i) - C7ARG; | ||
| AirMatrices.C7(i + 1) = AirMatrices.C7(i + 1) + C7ARG; | ||
| end | ||
| end |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,20 @@ | ||
| function [AirVariabes, RHS, SAVE, P_gg] = solveDryAirEquations(SoilVariables, GasDispersivity, TransportCoefficient, InitialValues, VaporVariables, ... | ||
| BoundaryCondition, ForcingData, P_gg, P_g, Xah, XaT, Xaa, RHODA, KT, Delt_t) | ||
| %{ | ||
| Solve the dry air equation with the Thomas algorithm to update the soil | ||
| air pressure 'P_gg', the finite difference time-stepping scheme is | ||
| exampled as for the soil moisture equation, which derived in 'STEMMUS | ||
| Technical Notes' section 4, Equation 4.32. | ||
| %} | ||
| AirVariabes = dryair.calculateDryAirParameters(SoilVariables, GasDispersivity, TransportCoefficient, InitialValues, VaporVariables, ... | ||
| P_gg, Xah, XaT, Xaa, RHODA); | ||
|
|
||
| AirMatrices = dryair.calculateMatricCoefficients(AirVariabes, InitialValues); | ||
|
|
||
| [RHS, AirMatrices, SAVE] = dryair.assembleCoefficientMatrices(AirMatrices, SoilVariables, Delt_t, P_g); | ||
|
|
||
| [RHS, AirMatrices] = dryair.calculateBoundaryConditions(BoundaryCondition, AirMatrices, ForcingData, RHS, KT); | ||
|
|
||
| [AirMatrices, P_gg, RHS] = dryair.solveTridiagonalMatrixEquations(RHS, AirMatrices); | ||
|
|
||
| end |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,22 @@ | ||
| function [AirMatrices, P_gg, RHS] = solveTridiagonalMatrixEquations(RHS, AirMatrices) | ||
| %{ | ||
| Use Thomas algorithm to solve the tridiagonal matrix equations, which is | ||
| in the form of Equation 4.25 STEMMUS Technical Notes, page 41. | ||
| %} | ||
|
|
||
| ModelSettings = io.getModelSettings(); | ||
| RHS(1) = RHS(1) / AirMatrices.C6(1, 1); | ||
|
|
||
| for i = 2:ModelSettings.NN | ||
| AirMatrices.C6(i, 1) = AirMatrices.C6(i, 1) - AirMatrices.C6_a(i - 1) * AirMatrices.C6(i - 1, 2) / AirMatrices.C6(i - 1, 1); | ||
| RHS(i) = (RHS(i) - AirMatrices.C6_a(i - 1) * RHS(i - 1)) / AirMatrices.C6(i, 1); | ||
| end | ||
|
|
||
| for i = ModelSettings.NL:-1:1 | ||
| RHS(i) = RHS(i) - AirMatrices.C6(i, 2) * RHS(i + 1) / AirMatrices.C6(i, 1); | ||
| end | ||
|
|
||
| for i = 1:ModelSettings.NN | ||
| P_gg(i) = RHS(i); | ||
| end | ||
| end |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Oops, something went wrong.