VBA Matrix Operations Help Needed.

[dennygibson]

Hi … here's a snippet from a macro I'm working on … I can't get the last statement (multiply a matrix by the inverse of another matrix) Everything else works but the last statement:


Option Explicit
Option Base 1
Sub Coefficients()
Dim C() As Variant
Dim D() As Variant
Dim F() As Variant
ReDim C(1 To 2, 1 To 2)
ReDim D(1 To 2, 1 To 1)
ReDim F(1 To 2, 1 To 1)
ReDim Cinv(1 To 2, 1 To 2)
'Set initial conditions
W1 = 0
Del1 = 0.05
' Swirl Chamber Coefficient Assembly
C(1, 1) = Del1 * W1
C(1, 2) = Del1 ^ 2
C(2, 1) = ((A2 - B2) * Del1 * W1 ^ 2)
C(2, 2) = ((A2 - 2 * B2 + 1) * Del1 ^ 2 * W1)
D(1, 1) = B1 / A1
D(2, 1) = C2 * W1
ActiveSheet.Range("T5") = C(1, 1)
ActiveSheet.Range("U5") = C(1, 2)
ActiveSheet.Range("T6") = C(2, 1)
ActiveSheet.Range("U6") = C(2, 2)
ActiveSheet.Range("X5") = D(1, 1)
ActiveSheet.Range("X6") = D(2, 1)
' Get F
' Note: we have {C}dx/dt={D} ... need: dx/dt={F}=[C_inverse]{D}
F = WorksheetFunction.MMult(WorksheetFunction.MInverse(C), D) --- this is the statement that bombs!!!
End Sub

[Paul_Hossler]

Welcome to the forum

Take a minute to read the links and the FAQ in my sigature, especiallt the part about using code tags and including a small sample workbook with the macros and some data that shows the issue or desired result


As a GUESS …




Option Explicit
Option Base 1
Sub Coefficients()
Dim C() As Variant
Dim D() As Variant
Dim F As Variant '--------------------------------------------- make F a Variant, not an array of Variants
ReDim C(1 To 2, 1 To 2)
ReDim D(1 To 2, 1 To 1)
' ReDim F(1 To 2, 1 To 1) ------------------------------------------------------------------- no need to ReDim since above change
ReDim Cinv(1 To 2, 1 To 2)
'Set initial conditions
W1 = 0
Del1 = 0.05
' Swirl Chamber Coefficient Assembly
C(1, 1) = Del1 * W1
C(1, 2) = Del1 ^ 2
C(2, 1) = ((A2 - B2) * Del1 * W1 ^ 2)
C(2, 2) = ((A2 - 2 * B2 + 1) * Del1 ^ 2 * W1)
D(1, 1) = B1 / A1
D(2, 1) = C2 * W1
ActiveSheet.Range("T5") = C(1, 1)
ActiveSheet.Range("U5") = C(1, 2)
ActiveSheet.Range("T6") = C(2, 1)
ActiveSheet.Range("U6") = C(2, 2)
ActiveSheet.Range("X5") = D(1, 1)
ActiveSheet.Range("X6") = D(2, 1)
' Get F
' Note: we have {C}dx/dt={D} ... need: dx/dt={F}=[C_inverse]{D}
F = WorksheetFunction.MMult(WorksheetFunction.MInverse(C), D) --- this is the statement that bombs!!!
End Sub

[MagPower]

How are A2, B2, C2, defined, and what are their values?

[dennygibson]

How are A2, B2, C2, defined, and what are their values?

I knew someone would ask … I didn't include them because the whole code is fairly long (and getting longer). Here's the hole thing so far:


Option Explicit
Option Base 1
Sub Coefficients()
Dim ShMain As Worksheet
Dim Rswrl As Variant
Dim Lswrl As Variant
Dim Ror As Variant
Dim Alpha As Variant
Dim rho_liq As Variant
Dim visc As Variant
Dim rho_air As Variant
Dim surfT_liquid As Variant
Dim Omega As Variant
Dim Re As Variant
Dim Ain As Variant
Dim V1 As Variant
Dim Ri As Variant
Dim R As Variant
Dim psi As Variant
Dim psipsi As Variant
Dim a As Variant
Dim b As Variant
Dim h As Variant
Dim sum As Variant
Dim First As Variant
Dim Last As Variant
Dim A1 As Variant
Dim A2 As Variant
Dim B1 As Variant
Dim B2 As Variant
Dim C2 As Variant
Dim A3 As Variant
Dim B3 As Variant
Dim B33 As Variant
Dim C3 As Variant
Dim B4 As Variant
Dim A4 As Variant
Dim W1 As Variant
Dim Del1 As Variant
Dim C As Variant
Dim D As Variant
Dim F As Variant
Dim Pi As Variant
ReDim C(1 To 2, 1 To 2)
ReDim D(1 To 2, 1 To 1)
ReDim F(1 To 2, 1 To 1)
ReDim Cinv(1 To 2, 1 To 2)
Dim npts As Long
Dim k As Long
' Input Values:
' Geometry:
Rswrl = ActiveSheet.Range("H18").Value / 2
Lswrl = ActiveSheet.Range("H14").Value
Ror = ActiveSheet.Range("H20").Value / 2
Alpha = ActiveSheet.Range("H21").Value
Ri = ActiveSheet.Range("H17").Value
' Properties:
rho_liq = ActiveSheet.Range("H8").Value
visc = ActiveSheet.Range("H9").Value
rho_air = ActiveSheet.Range("H10").Value
surfT_liquid = ActiveSheet.Range("H11").Value
' Preliminaries
Pi = 3.14159265
Ain = ActiveSheet.Range("H23").Value
V1 = ActiveSheet.Range("H24").Value
' Swirl Chamber
Omega = V1 * (Rswrl - Ri)
Re = Omega / visc
R = Rswrl
' Coefficients (constants - variable values must be added.
' Coefficient A1
npts = 100
b = 1#
a = 0#
h = (b - a) / npts

' First and Last Integration Points
First = (2 * a - a ^ 2) - (2 * a - a ^ 2) ^ 2
Last = (2 * b - b ^ 2) - (2 * b - b ^ 2) ^ 2
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 4 * (psi - psi * psi)
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 2 * (psi - psi * psi)
Next k
' Add final value
sum = sum + Last
A1 = h / 3 * sum
ActiveSheet.Range("R5") = A1

' Coefficient A2
npts = 50
b = 1#
a = 0#
h = (b - a) / npts
' First and Last Integration Points
First = (2 * a - a ^ 2)
Last = (2 * b - b ^ 2)
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 4 * psi
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 2 * psi
Next k
' Add final value
sum = sum + Last
A2 = h / 3 * sum
ActiveSheet.Range("R6") = A2

' Coefficient B1 and C2
b = 1#
a = 0#
B1 = (2 - 2 * b) - (2 - 2 * a)
ActiveSheet.Range("R7") = B1
C2 = B1
ActiveSheet.Range("R9") = C2

' Coefficient B2
npts = 100
b = 1#
a = 0#
h = (b - a) / npts
' First and Last Integration Points
First = (2 * a - a ^ 2) ^ 2
Last = (2 * b - b ^ 2) ^ 2
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psipsi = (2 * k * h - (k * h) * (k * h)) * (2 * k * h - (k * h) * (k * h))
sum = sum + 4 * psipsi
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psipsi = (2 * k * h - (k * h) * (k * h)) * (2 * k * h - (k * h) * (k * h))
sum = sum + 2 * psipsi
Next k
' Add final value
sum = sum + Last
B2 = h / 3 * sum
ActiveSheet.Range("R8") = B2
'Set initial conditions
W1 = 0
Del1 = 0.05
' Swirl Chamber Coefficient Assembly
C(1, 1) = Del1 * W1
C(1, 2) = Del1 ^ 2
C(2, 1) = ((A2 - B2) * Del1 * W1 ^ 2)
C(2, 2) = ((A2 - 2 * B2 + 1) * Del1 ^ 2 * W1)
D(1, 1) = B1 / A1
D(2, 1) = C2 * W1
ActiveSheet.Range("T5") = C(1, 1)
ActiveSheet.Range("U5") = C(1, 2)
ActiveSheet.Range("T6") = C(2, 1)
ActiveSheet.Range("U6") = C(2, 2)
ActiveSheet.Range("X5") = D(1, 1)
ActiveSheet.Range("X6") = D(2, 1)
' Get F
F = WorksheetFunction.MMult(WorksheetFunction.MInverse(C), D)


' Spin Chamber
' Coefficients
A3 = A1
B3 = A2
B33 = B2
C3 = B1
B4 = C2
A4 = A1
ActiveSheet.Range("R10") = A3
ActiveSheet.Range("R11") = B3
ActiveSheet.Range("R12") = B33
ActiveSheet.Range("R13") = C3
ActiveSheet.Range("R14") = B4
ActiveSheet.Range("R15") = A4

' Spray Orifice
Omega = ActiveSheet.Range("H25")
Re = ActiveSheet.Range("H26")
R = Ror

' Solve Momentum Equations


' Zone 1

End Sub

[dennygibson]

Paul … thanks! I read it. Ok … I'm not familiar with [/CODE] tags rather than ask you what they do I'll go add to my code and see for myself. Also … I'm a bit confused … F is a 2x1 matrix.

[Paul_Hossler]

1. You can use the [#] icon to insert CODE and /CODE tags -- paste the macro between then

2. The reason I suggested attaching a workbook is you're getting values from the worksheet



Rswrl = ActiveSheet.Range("H18").Value / 2
Lswrl = ActiveSheet.Range("H14").Value
Ror = ActiveSheet.Range("H20").Value / 2
Alpha = ActiveSheet.Range("H21").Value
Ri = ActiveSheet.Range("H17").Value

' Properties:
rho_liq = ActiveSheet.Range("H8").Value
visc = ActiveSheet.Range("H9").Value
rho_air = ActiveSheet.Range("H10").Value
surfT_liquid = ActiveSheet.Range("H11").Value


and I don't have any idea what a reasonable guess-timent is. Bottom right, [Go Advanced] and use the paperclip icon

3. When I ran the macro with my guesses, I go the error on MInverse, and a Watch window shows.

23740

I'm fuzzy, but I don't think C is invertable (key word is 'fuzzy')

[dennygibson]

Paul … when I pasted my code, highlighted it, and pressed # … the code disappeared and left only # in its place. Also … your suggestion didn't fix it. Here's my code again with better formatting:

Dim visc As Variant
Dim rho_air As Variant
Dim surfT_liquid As Variant
Dim Omega As Variant
Dim Re As Variant
Dim Ain As Variant
Dim V1 As Variant
Dim Ri As Variant
Dim R As Variant
Dim psi As Variant
Dim psipsi As Variant
Dim a As Variant
Dim b As Variant
Dim h As Variant
Dim sum As Variant
Dim First As Variant
Dim Last As Variant
Dim A1 As Variant
Dim A2 As Variant
Dim B1 As Variant
Dim B2 As Variant
Dim C2 As Variant
Dim A3 As Variant
Dim B3 As Variant
Dim B33 As Variant
Dim C3 As Variant
Dim B4 As Variant
Dim A4 As Variant
Dim W1 As Variant
Dim Del1 As Variant
Dim C()As Variant
Dim D() As Variant
Dim F As Variant
Dim Pi As Variant
ReDim C(1 To 2, 1 To 2)
ReDim D(1 To 2, 1 To 1)
'ReDim F(1 To 2, 1 To 1)
Dim npts As Long
Dim k As Long
' Input Values:
' Geometry:
Rswrl = ActiveSheet.Range("H18").Value / 2
Lswrl = ActiveSheet.Range("H14").Value
Ror = ActiveSheet.Range("H20").Value / 2
Alpha = ActiveSheet.Range("H21").Value
Ri = ActiveSheet.Range("H17").Value
' Properties:
rho_liq = ActiveSheet.Range("H8").Value
visc = ActiveSheet.Range("H9").Value
rho_air = ActiveSheet.Range("H10").Value
surfT_liquid = ActiveSheet.Range("H11").Value
' Preliminaries
Pi = 3.14159265
Ain = ActiveSheet.Range("H23").Value
V1 = ActiveSheet.Range("H24").Value
' Swirl Chamber
Omega = V1 * (Rswrl - Ri)
Re = Omega / visc
R = Rswrl
' Coefficients (constants - variable values must be added.
' Coefficient A1
npts = 100
b = 1#
a = 0#
h = (b - a) / npts

' First and Last Integration Points
First = (2 * a - a ^ 2) - (2 * a - a ^ 2) ^ 2
Last = (2 * b - b ^ 2) - (2 * b - b ^ 2) ^ 2
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 4 * (psi - psi * psi)
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 2 * (psi - psi * psi)
Next k
' Add final value
sum = sum + Last
A1 = h / 3 * sum
ActiveSheet.Range("R5") = A1

' Coefficient A2
npts = 50
b = 1#
a = 0#
h = (b - a) / npts
' First and Last Integration Points
First = (2 * a - a ^ 2)
Last = (2 * b - b ^ 2)
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 4 * psi
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 2 * psi
Next k
' Add final value
sum = sum + Last
A2 = h / 3 * sum
ActiveSheet.Range("R6") = A2

' Coefficient B1 and C2
b = 1#
a = 0#
B1 = (2 - 2 * b) - (2 - 2 * a)
ActiveSheet.Range("R7") = B1
C2 = B1
ActiveSheet.Range("R9") = C2

' Coefficient B2
npts = 100
b = 1#
a = 0#
h = (b - a) / npts
' First and Last Integration Points
First = (2 * a - a ^ 2) ^ 2
Last = (2 * b - b ^ 2) ^ 2
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psipsi = (2 * k * h - (k * h) * (k * h)) * (2 * k * h - (k * h) * (k * h))
sum = sum + 4 * psipsi
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psipsi = (2 * k * h - (k * h) * (k * h)) * (2 * k * h - (k * h) * (k * h))
sum = sum + 2 * psipsi
Next k
' Add final value
sum = sum + Last
B2 = h / 3 * sum
ActiveSheet.Range("R8") = B2
'Set initial conditions
W1 = 0
Del1 = 0.05
' Swirl Chamber Coefficient Assembly
C(1, 1) = Del1 * W1
C(1, 2) = Del1 ^ 2
C(2, 1) = ((A2 - B2) * Del1 * W1 ^ 2)
C(2, 2) = ((A2 - 2 * B2 + 1) * Del1 ^ 2 * W1)
D(1, 1) = B1 / A1
D(2, 1) = C2 * W1
ActiveSheet.Range("T5") = C(1, 1)
ActiveSheet.Range("U5") = C(1, 2)
ActiveSheet.Range("T6") = C(2, 1)
ActiveSheet.Range("U6") = C(2, 2)
ActiveSheet.Range("X5") = D(1, 1)
ActiveSheet.Range("X6") = D(2, 1)
' Get F
F = WorksheetFunction.MMult(WorksheetFunction.MInverse(C), D) 'this is the problem statement.


' Spin Chamber
' Coefficients
A3 = A1
B3 = A2
B33 = B2
C3 = B1
B4 = C2
A4 = A1
ActiveSheet.Range("R10") = A3
ActiveSheet.Range("R11") = B3
ActiveSheet.Range("R12") = B33
ActiveSheet.Range("R13") = C3
ActiveSheet.Range("R14") = B4
ActiveSheet.Range("R15") = A4

' Spray Orifice
Omega = ActiveSheet.Range("H25")
Re = ActiveSheet.Range("H26")
R = Ror

' Solve Momentum Equations


' Zone 1

End Sub

[MagPower]

Just curious, is this a problem from fluid mechanics? :think:

[dennygibson]

Paul … yes … it's the beginnings of a theoretical, viscous simulation of a pressure swirl atomizer. Let me see if I can get this right this time ... I attached the workbook … hopefully you can access it. I changed the initial conditions to small numbers so that each cell in the matrix had a value and got the inverse to write to the worksheet. Apparently you can't invert a matrix with zeros in all but one array location (guess it makes sense).

Now I'm on to figuring out Runge Kutta 4th order with matrices … not a clue yet.

[dennygibson]

and let me try to post the code with tags:


Option Explicit
Option Base 1
Sub Coefficients()
Dim Rswrl As Variant
Dim Lswrl As Variant
Dim Ror As Variant
Dim Alpha As Variant
Dim rho_liq As Variant
Dim visc As Variant
Dim rho_air As Variant
Dim surfT_liquid As Variant
Dim Omega As Variant
Dim Re As Variant
Dim Ain As Variant
Dim V1 As Variant
Dim Ri As Variant
Dim R As Variant
Dim psi As Variant
Dim psipsi As Variant
Dim a As Variant
Dim b As Variant
Dim h As Variant
Dim sum As Variant
Dim First As Variant
Dim Last As Variant
Dim A1 As Variant
Dim A2 As Variant
Dim B1 As Variant
Dim B2 As Variant
Dim C2 As Variant
Dim A3 As Variant
Dim B3 As Variant
Dim B33 As Variant
Dim C3 As Variant
Dim B4 As Variant
Dim A4 As Variant
Dim W1 As Variant
Dim Del1 As Variant
Dim C()As Variant
Dim D() As Variant
Dim F As Variant
Dim Pi As Variant
ReDim C(1 To 2, 1 To 2)
ReDim D(1 To 2, 1 To 1)
'ReDim F(1 To 2, 1 To 1)
Dim npts As Long
Dim k As Long
' Input Values:
' Geometry:
Rswrl = ActiveSheet.Range("H18").Value / 2
Lswrl = ActiveSheet.Range("H14").Value
Ror = ActiveSheet.Range("H20").Value / 2
Alpha = ActiveSheet.Range("H21").Value
Ri = ActiveSheet.Range("H17").Value
' Properties:
rho_liq = ActiveSheet.Range("H8").Value
visc = ActiveSheet.Range("H9").Value
rho_air = ActiveSheet.Range("H10").Value
surfT_liquid = ActiveSheet.Range("H11").Value
' Preliminaries
Pi = 3.14159265
Ain = ActiveSheet.Range("H23").Value
V1 = ActiveSheet.Range("H24").Value
' Swirl Chamber
Omega = V1 * (Rswrl - Ri)
Re = Omega / visc
R = Rswrl
' Coefficients (constants - variable values must be added.
' Coefficient A1
npts = 100
b = 1#
a = 0#
h = (b - a) / npts

' First and Last Integration Points
First = (2 * a - a ^ 2) - (2 * a - a ^ 2) ^ 2
Last = (2 * b - b ^ 2) - (2 * b - b ^ 2) ^ 2
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 4 * (psi - psi * psi)
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 2 * (psi - psi * psi)
Next k
' Add final value
sum = sum + Last
A1 = h / 3 * sum
ActiveSheet.Range("R5") = A1

' Coefficient A2
npts = 50
b = 1#
a = 0#
h = (b - a) / npts
' First and Last Integration Points
First = (2 * a - a ^ 2)
Last = (2 * b - b ^ 2)
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 4 * psi
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 2 * psi
Next k
' Add final value
sum = sum + Last
A2 = h / 3 * sum
ActiveSheet.Range("R6") = A2

' Coefficient B1 and C2
b = 1#
a = 0#
B1 = (2 - 2 * b) - (2 - 2 * a)
ActiveSheet.Range("R7") = B1
C2 = B1
ActiveSheet.Range("R9") = C2

' Coefficient B2
npts = 100
b = 1#
a = 0#
h = (b - a) / npts
' First and Last Integration Points
First = (2 * a - a ^ 2) ^ 2
Last = (2 * b - b ^ 2) ^ 2
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psipsi = (2 * k * h - (k * h) * (k * h)) * (2 * k * h - (k * h) * (k * h))
sum = sum + 4 * psipsi
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psipsi = (2 * k * h - (k * h) * (k * h)) * (2 * k * h - (k * h) * (k * h))
sum = sum + 2 * psipsi
Next k
' Add final value
sum = sum + Last
B2 = h / 3 * sum
ActiveSheet.Range("R8") = B2
'Set initial conditions
W1 = 0
Del1 = 0.05
' Swirl Chamber Coefficient Assembly
C(1, 1) = Del1 * W1
C(1, 2) = Del1 ^ 2
C(2, 1) = ((A2 - B2) * Del1 * W1 ^ 2)
C(2, 2) = ((A2 - 2 * B2 + 1) * Del1 ^ 2 * W1)
D(1, 1) = B1 / A1
D(2, 1) = C2 * W1
ActiveSheet.Range("T5") = C(1, 1)
ActiveSheet.Range("U5") = C(1, 2)
ActiveSheet.Range("T6") = C(2, 1)
ActiveSheet.Range("U6") = C(2, 2)
ActiveSheet.Range("X5") = D(1, 1)
ActiveSheet.Range("X6") = D(2, 1)
' Get F
F = WorksheetFunction.MMult(WorksheetFunction.MInverse(C), D)


' Spin Chamber
' Coefficients
A3 = A1
B3 = A2
B33 = B2
C3 = B1
B4 = C2
A4 = A1
ActiveSheet.Range("R10") = A3
ActiveSheet.Range("R11") = B3
ActiveSheet.Range("R12") = B33
ActiveSheet.Range("R13") = C3
ActiveSheet.Range("R14") = B4
ActiveSheet.Range("R15") = A4

' Spray Orifice
Omega = ActiveSheet.Range("H25")
Re = ActiveSheet.Range("H26")
R = Ror

' Solve Momentum Equations


' Zone 1

End Sub

[dennygibson]

MagPower ... sorry ... I got your reply mixed up with Paul's ... yes, this is a fluid mechnics problem ... trying to simulate swirling flow in a simplex atomizer.

[Paul_Hossler]

@dennygibson --

+1 on the CODE tags :thumb

+1 on attaching WB with 'real' data:thumb:thumb. I didn't understand it, but I'll just go with the flow :)

+1 on running to completion with no debug messages

Few minor suggestions about declaring variables you might consider

a. Since they're numbers, a Double is more efficient that a Variant
b. Sum is a Excel keyword, so a different one might be less confusing
c .You can Dim arrays with dimensions if the size isn't going to change
d. You can use Excel's Pi




Dim Rswrl As Double ' Double is more efficient than Variants if numbers
Dim Lswrl As Double
Dim Ror As Double
Dim Alpha As Double
Dim rho_liq As Double
Dim visc As Double
Dim rho_air As Double
Dim surfT_liquid As Double
Dim Omega As Double
Dim Re As Double
Dim Ain As Double
Dim V1 As Double
Dim Ri As Double
Dim R As Double
Dim psi As Double
Dim psipsi As Double
Dim a As Double
Dim b As Double
Dim h As Double
Dim sum As Double ' Sum is keyword - different one would be better
Dim First As Double
Dim Last As Double
Dim A1 As Double
Dim A2 As Double
Dim B1 As Double
Dim B2 As Double
Dim C2 As Double
Dim A3 As Double
Dim B3 As Double
Dim B33 As Double
Dim C3 As Double
Dim B4 As Double
Dim A4 As Double
Dim W1 As Double
Dim Del1 As Double
Dim C(1 To 2, 1 To 2) As Double ' no need to ReDim - just size at beginning
Dim D(1 To 2, 1 To 1) As Double
Dim C_inv As Variant ' these have to be Variant since you assign to them
Dim F As Variant
Dim Pi As Double
Pi = Application.WorksheetFunction.Pi ' you can use the built in or make a Const
'Const Pi as double = 3.14159265358979


But what the heck … it works :yes

[dennygibson]

Paul ... thanks! I appreciate the suggestions. Would this be the place to get some help using runge-kutta 4th order to solve simultaneous equations (using matrix operations) or am I barking up the wrong tree?

[MagPower]

Paul ... thanks! I appreciate the suggestions. Would this be the place to get some help using runge-kutta 4th order to solve simultaneous equations (using matrix operations) or am I barking up the wrong tree?

Denny,

I don't know if anyone in this forum has a background using RK to solve a system of DEs. A math forum might be a better place.

-- Russ

[dennygibson]

Russ ... thanks. I’ll see what I can find.

[MagPower]

Russ ... thanks. I’ll see what I can find.

Denny - I found this site. Hopefully this helps:

https://www.researchgate.net/publication/259079149_A_Spreadsheet_Solution_of_a_System_of_Ordinary_Differential_Equat ions_Using_the_Fourth-Order_Runge-Kutta_Method

Russ

[dennygibson]

Russ ... thanks! I’ll check it out

[dennygibson]

Hey! Things are going well but I have a variable type mismatch that I haven't been able to figure out. The three statements are highlighted in red. I'm not done compiling but I got hung up there. The workbook is attached here 23752 Take a look and save me again?



Option Explicit
Option Base 1
Sub Coefficients()
Dim Rswrl As Double
Dim Lswrl As Double
Dim Ror As Double
Dim Alpha As Double
Dim rho_liq As Double
Dim visc As Double
Dim rho_air As Double
Dim surfT_liquid As Double
Dim Omega As Double
Dim Re As Double
Dim Ain As Double
Dim V1 As Double
Dim Ri As Double
Dim R As Double
Dim psi As Double
Dim psipsi As Double
Dim a As Double
Dim b As Double
Dim h As Double
Dim sum As Double
Dim First As Double
Dim Last As Double
Dim A1 As Double
Dim A2 As Double
Dim B1 As Double
Dim B2 As Double
Dim C2 As Double
Dim A3 As Double
Dim B3 As Double
Dim B33 As Double
Dim C3 As Double
Dim B4 As Double
Dim A4 As Double
Dim C(1 To 2, 1 To 2) As Double
Dim D(1 To 2, 1 To 1) As Variant
Dim C_inv As Variant
Dim F As Variant
Dim Zswrl As Double
Dim dz1 As Double
Dim K1 As Double
Dim K2 As Double
Dim K3 As Double
Dim K4 As Double
Dim L1 As Double
Dim L2 As Double
Dim L3 As Double
Dim L4 As Double
Dim F1 As Double
Dim F2 As Double
Dim Z1 As Double
Dim Z2 As Double
Dim Z3 As Double
Dim Z4 As Double
Dim Del1 As Variant
Dim Del2 As Double
Dim Del3 As Double
Dim Del4 As Double
Dim W1 As Variant
Dim W2 As Double
Dim W3 As Double
Dim W4 As Double
Dim Del1St As Variant
Dim W1St As Variant
Dim Pi As Double
Dim npts As Long
Dim k As Long
Dim Nswrl As Long
Dim IZ As Long
' Input Values:
' Geometry:
Rswrl = ActiveSheet.Range("H18").Value / 2
Lswrl = ActiveSheet.Range("H14").Value
Ror = ActiveSheet.Range("H20").Value / 2
Alpha = ActiveSheet.Range("H21").Value
Ri = ActiveSheet.Range("H17").Value
' Properties:
rho_liq = ActiveSheet.Range("H8").Value
visc = ActiveSheet.Range("H9").Value
rho_air = ActiveSheet.Range("H10").Value
surfT_liquid = ActiveSheet.Range("H11").Value
' Preliminaries
Pi = 3.14159265
Ain = ActiveSheet.Range("H23").Value
V1 = ActiveSheet.Range("H24").Value
' Swirl Chamber
' Coefficients (constants - variable values must be added.
' Coefficient A1
npts = 20
b = 1#
a = 0#
h = (b - a) / npts

' First and Last Integration Points
First = (2 * a - a ^ 2) - (2 * a - a ^ 2) ^ 2
Last = (2 * b - b ^ 2) - (2 * b - b ^ 2) ^ 2
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 4 * (psi - psi * psi)
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 2 * (psi - psi * psi)
Next k
' Add final value
sum = sum + Last
A1 = h / 3 * sum
ActiveSheet.Range("R5") = A1

' Coefficient A2
npts = 20
b = 1#
a = 0#
h = (b - a) / npts
' First and Last Integration Points
First = (2 * a - a ^ 2)
Last = (2 * b - b ^ 2)
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 4 * psi
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psi = 2 * k * h - (k * h) * (k * h)
sum = sum + 2 * psi
Next k
' Add final value
sum = sum + Last
A2 = h / 3 * sum
ActiveSheet.Range("R6") = A2

' Coefficient B1 and C2
b = 1#
a = 0#
B1 = (2 - 2 * b) - (2 - 2 * a)
ActiveSheet.Range("R7") = B1
C2 = B1
ActiveSheet.Range("R9") = C2

' Coefficient B2
npts = 20
b = 1#
a = 0#
h = (b - a) / npts
' First and Last Integration Points
First = (2 * a - a ^ 2) ^ 2
Last = (2 * b - b ^ 2) ^ 2
sum = First
' Odd terms
For k = 1 To npts - 1 Step 2
psipsi = (2 * k * h - (k * h) * (k * h)) * (2 * k * h - (k * h) * (k * h))
sum = sum + 4 * psipsi
Next k
' Even terms
For k = 2 To npts - 1 Step 2
psipsi = (2 * k * h - (k * h) * (k * h)) * (2 * k * h - (k * h) * (k * h))
sum = sum + 2 * psipsi
Next k
' Add final value
sum = sum + Last
B2 = h / 3 * sum
ActiveSheet.Range("R8") = B2

' Zone 1
' Circulation Constant
Omega = V1 * (Rswrl - Ri)
' Reynolds Number
Re = Omega / visc
' Initial value
W1 = 0.001
Del1 = 0.001
Z1 = 0
' Swirl Chamber Height
Zswrl = ActiveSheet.Range("H18").Value
' Number of Steps in Swirl Chamber Height
Nswrl = 20
' Delta-Z
dz1 = Zswrl / Nswrl
For IZ = 1 To Nswrl
' Swirl Chamber Coefficient Matrix
C(1, 1) = Del1 * W1
C(1, 2) = Del1 ^ 2
C(2, 1) = ((A2 - B2) * Del1 * W1 ^ 2)
C(2, 2) = ((A2 - 2 * B2 + 1) * Del1 ^ 2 * W1)
D(1, 1) = B1 / A1
D(2, 1) = C2 * W1

' Vector F Initial Conditions
F = WorksheetFunction.MMult(WorksheetFunction.MInverse(C), D)
ActiveSheet.Range("Y5") = F(1, 1)
ActiveSheet.Range("Y6") = F(2, 1)
F1 = F(1, 1)
F2 = F(2, 1)
K1 = dz1 * F1
L1 = dz1 * F2
Z1 = Z1 + dz1 / 2
Del1 = Del1 + K1 / 2
W1 = W1 + L1 / 2
C(1, 1) = Del1 * W1
C(1, 2) = Del1 ^ 2
C(2, 1) = ((A2 - B2) * Del1 * W1 ^ 2)
C(2, 2) = ((A2 - 2 * B2 + 1) * Del1 ^ 2 * W1)
D(1, 1) = B1 / A1
D(2, 1) = C2 * W1
F = WorksheetFunction.MMult(WorksheetFunction.MInverse(C), D)
F1 = F(1, 1)
F2 = F(2, 1)
K2 = dz1 * F1
L2 = dz1 * F2
Z2 = Z1 + dz1 / 2
Del2 = Del1 + K2 / 2
W2 = W1 + L2 / 2
C(1, 1) = Del2 * W2
C(1, 2) = Del2 ^ 2
C(2, 1) = ((A2 - B2) * Del2 * W2 ^ 2)
C(2, 2) = ((A2 - 2 * B2 + 1) * Del2 ^ 2 * W2)
D(1, 1) = B1 / A1
D(2, 1) = C2 * W2
F = WorksheetFunction.MMult(WorksheetFunction.MInverse(C), D)
F1 = F(1, 1)
F2 = F(2, 1)
K3 = dz1 * F1
L3 = dz1 * F2
Z3 = Z1 + dz1 / 2
Del3 = Del1 + K2 / 2
W3 = W1 + L2 / 2
C(1, 1) = Del3 * W3
C(1, 2) = Del3 ^ 2
C(2, 1) = ((A2 - B2) * Del3 * W3 ^ 2)
C(2, 2) = ((A2 - 2 * B2 + 1) * Del3 ^ 2 * W3)
D(1, 1) = B1 / A1
D(2, 1) = C2 * W3
F = WorksheetFunction.MMult(WorksheetFunction.MInverse(C), D)
F1 = F(1, 1)
F2 = F(2, 1)
K4 = dz1 * F1
L4 = dz1 * F2
Z3 = Z1 + dz1 / 2
Del4 = Del1 + K3
W4 = W1 + L3
C(1, 1) = Del4 * W4
C(1, 2) = Del4 ^ 2
C(2, 1) = ((A2 - B2) * Del4 * W4 ^ 2)
C(2, 2) = ((A2 - 2 * B2 + 1) * Del4 ^ 2 * W4)
D(1, 1) = B1 / A1
D(2, 1) = C2 * W4

Del1 = Del1 + (K1 + 2 * K2 + 2 * K3 + K4) / 6
W1 = W1 + (L1 + 2 * L2 + 2 * L3 + L4) / 6
Del1St(IZ) = Del1 'Type Mismatch
W1St(IZ) = W1 'Type Mismatch
Z1St(IZ) = Z1 'Type Mismatch
' Adjust Z1
Z1 = Z1 + dz1
Next IZ

' Spin Chamber
' Coefficients
A3 = A1
B3 = A2
B33 = B2
C3 = B1
B4 = C2
A4 = A1
ActiveSheet.Range("R10") = A3
ActiveSheet.Range("R11") = B3
ActiveSheet.Range("R12") = B33
ActiveSheet.Range("R13") = C3
ActiveSheet.Range("R14") = B4
ActiveSheet.Range("R15") = A4

' Spray Orifice
Omega = ActiveSheet.Range("H25")
Re = ActiveSheet.Range("H26")
R = Ror

' Solve Momentum Equations


End Sub

[Paul_Hossler]

You Dim-ed Del1St and W1St as Variant, and plain variants don't take an Index (IZ)

I did not see where Z1St was Dim-ed at all





Del1St(IZ) = Del1 'Type Mismatch
W1St(IZ) = W1 'Type Mismatch
Z1St(IZ) = Z1 'Type Mismatch


I'm guessing you wanted







Dim Del1St(1 to 20) As Double
Dim W11St(1 to 20) As Double
Dim Z1St(1 to 20) As Double

[dennygibson]

Paul ... thanks! As usual you were very helpful!

[dennygibson]

Is there a way for the array range to be dynamic … like Dim array(1 to N) Double N could be set differently in different parts of the code?

[dennygibson]

Is there a way for the array range to be dynamic … like Dim array(1 to N) Double where N could be set differently in different parts of the code?

[dennygibson]

You Dim-ed Del1St and W1St as Variant, and plain variants don't take an Index (IZ)

I did not see where Z1St was Dim-ed at all





Del1St(IZ) = Del1 'Type Mismatch
W1St(IZ) = W1 'Type Mismatch
Z1St(IZ) = Z1 'Type Mismatch


I'm guessing you wanted






Dim Del1St(1 to 20) As Double
Dim W11St(1 to 20) As Double
Dim Z1St(1 to 20) As Double

[dennygibson]

Paul ... please help me remember how to mark a topic “solved”

[Paul_Hossler]

Above your first post, there's [Thread Tools] and an option is to mark [SOLVED]

[Paul_Hossler]

Is there a way for the array range to be dynamic … like Dim array(1 to N) Double where N could be set differently in different parts of the code?


Yes

Look in Help for "ReDim" and "Redim Preserve"

However, I've found that most times it's not needed and can be designed around

Some sample code




Option Explicit

Sub Arrays()
Dim A(1 To 2) As Long ' fixed
Dim B() As Long ' ReDim-able

A(1) = 100
A(2) = 200
MsgBox A(1) & " -- " & A(2)

ReDim B(1 To 3)
B(1) = 1000
B(2) = 2000
B(3) = 3000
MsgBox B(1) & " -- " & B(2) & " -- " & B(3)

ReDim Preserve B(1 To 4)
B(4) = 4000
MsgBox B(1) & " -- " & B(2) & " -- " & B(3) & " -- " & B(4)

Erase A
MsgBox A(1) & " -- " & A(2)
A(1) = 10000
A(2) = 20000
MsgBox A(1) & " -- " & A(2)
End Sub

[dennygibson]

Paul ... thanks!