How to solve a Matrix of ODEs

Question

Hey could you help me solving this matrix of linear ODEs. This exercise comes from the book "Transport Phenomena" by Byrd, Stewart and Lightfoot [23.D.4]. I managed to get this matrix of ODEs, i tried solving it by hand but i can't get to the solution.

$$ \frac{d(C_2)}{dt} =h\cdot C_1-B\cdot C_2 $$

Where:

$$ C_0 = \begin{bmatrix} C_{A0} \\ C_{B0} \\ C_{C0} \\ \end{bmatrix} \quad \:\:\: ; \:\:\:\:\:\:\: C_1 = \begin{bmatrix} C_{A1} \\ C_{B1} \\ C_{C1} \\ \end{bmatrix} \quad \:\:\: ; \:\:\:\:\:\:\: C_2 = \begin{bmatrix} C_{A2} \\ C_{B2} \\ C_{C2} \\ \end{bmatrix} \quad $$

And

$$ B= \begin{bmatrix} (h+k_{1}) & -k_{1} & 0 \\ - k_{1} & (h+k'_{1}+k_{2}) & -k_{2}\\ 0 & -k_{2} & (h+k'_{2})\\ \end{bmatrix} $$

Also, by solving other ODE I managed to get this result, and I know this is correct:

$$ C_1 =h(I-e^{-Bt})B^{-1}\cdot C_{0} $$

According to the book the solution is: $$ C_2 =h^2(I-e^{-Bt}-Bte^{-Bt})B^{-2}\cdot C_{0} $$

Saddly I've tried different "handmade" attempts but i can't get this result. This is what I always end up with: $$ C_2 =h^2(I-Bte^{-Bt})B^{-2}\cdot C_{0} $$

It seems to me that there is no use in writing the values of the B matrix in the code since in the result it only appears as B. I also tried using an integration factor since it is a linear ODE, but i just keep getting the same result as earlier

This is what I wrote in Mathematica, but i could't get anywhere (Im not good at writing code, don't hate me):

B = {{k1 + h, -k1, 0}, {-k1, 
h + k1´ + k2, -k2´}, {0, -k2, h + k2´}} // MatrixForm;
C0 = {{CA0}, {CB0}, {CC0}} // MatrixForm;
C1[t_] = h.(IdentityMatrix - Exp[-B*t])*(B^-1)*C0 // MatrixForm;
C2[t_] = {CA2[t], CB2[t], CC2[t]} // MatrixForm;
sol = DSolve[C2'[t] == h.C1[t] - B.C2'[t], C2[t], t]

Compare this, character by character, with what you have and see if this fixes the first couple of layers of misunderstandings. B={{k1+h,-k1,0},{-k1,h+k1+k2,-k2},{0,-k2,h+k2}}; C0={CA0,CB0,CC0}; C1[t_]:=h.(IdentityMatrix[3]-Exp[-B*t]).Inverse[B].C0; C2[t_]:={CA2[t],CB2[t],CC2[t]}; sol=DSolve[C2'[t]==h.C1[t]-B.C2'[t],C2[t],t] More will be needed but perhaps this is a start. See if you can use the documentation to understand why each change might have been made. It worries me that you defined the function C2[t] and then you want DSolve to solve for C2[t] — Bill, Commented Jul 29, 2020 at 19:08
Thanks that actually makes more sense. I still get a result that is way too long to be the correct result, but i'll keep trying. Maybe im wrong at how im using Dsolve, i dont really now how to use it. — felipe, Commented Jul 29, 2020 at 19:15
It would be a good idea to remove // MatrixForm from everything - it's just a presentational form — flinty, Commented Jul 29, 2020 at 19:15
Maybe what you want is DSolve[C2'[t]==h.C1[t]-B.C2'[t],{CA2[t], CB2[t], CC2[t]},t] — Bill, Commented Jul 29, 2020 at 19:17
What are k1,k2 and k1',k2' supposed to mean? Are the ' derivatives and are the k1,k2 functions of time? — flinty, Commented Jul 29, 2020 at 21:46

Steffen Jaeschke · Accepted Answer · 2020-07-29 21:48:38Z

This seems to be a hasty start into Mathematica. Your inputs are right but used without insight.

Use

B = {{k1 + h, -k1, 0}, {-k1,h + k1´ + k2, -k2´}, {0, -k2, h + k2´}}; C0 = {{CA0}, {CB0}, {CC0}}; C1[t_] = h.(IdentityMatrix - Exp[-Bt])(B^-1)*C0; C2[t_] = {CA2[t], CB2[t], CC2[t]};

B// MatrixForm

C0 // MatrixForm

C1[t] // MatrixForm

C2[t] // MatrixForm

This is the differential equation.

{C1[0],C2[0]}

is another check. This shows Mathematica is not omniscient towards user input.

What is dividing by a matrix in mathematics?

It is multiplying with the Inverse matrix. Now Mathematica can calculate the arbitrary inverse of symbolic matrices but for the given one it does not.

So the question is transformed into one more precise one. Help is required in how to enter a matrix for which this is calculated in Mathematica.

Most probable is that Mathematica expects the k1´, k2´ to be derivatives. So the suggestion of Bill is pointing in the right direction. The ´ in this case denotes a different set of reaction constant coefficients in a two-step flow reaction. The ´ can not be neglected because they are different. So simply rename them.

B = {{k1 + h, -k1, 0}, {-k1, h + kp1 + k2, -kp2}, {0, -k2, h + kp2}}

The next point in discourse is that naming the constants k in Mathematica leads to confusion because Wolfram Inc did choose k for naming free constants or even integration variables.

Keeping in simpler leave aside the meaning of the matrix coefficients. Another bigger problem is that Mathematica does not give such nice formulas as the introductory book does as solution and that is really worse in this case.

Det[B]

(-k1^2 + (h + k1) (h + k1\.b4 + k2)) (h + k2\.b4) + 
 k2 (-h k2\.b4 - k1 k2\.b4)

Inverse[B]

{{(h^2 + h k1\.b4 + h k2 + h k2\.b4 + 
   k1\.b4 k2\.b4)/((-k1^2 + (h + k1) (h + k1\.b4 + k2)) (h + k2\.b4) +
    k2 (-h k2\.b4 - k1 k2\.b4)), (
  h k1 + k1 k2\.b4)/((-k1^2 + (h + k1) (h + k1\.b4 + k2)) (h + 
      k2\.b4) + k2 (-h k2\.b4 - k1 k2\.b4)), (
  k1 k2\.b4)/((-k1^2 + (h + k1) (h + k1\.b4 + k2)) (h + k2\.b4) + 
   k2 (-h k2\.b4 - k1 k2\.b4))}, {(
  h k1 + k1 k2\.b4)/((-k1^2 + (h + k1) (h + k1\.b4 + k2)) (h + 
      k2\.b4) + 
   k2 (-h k2\.b4 - k1 k2\.b4)), ((h + k1) (h + 
     k2\.b4))/((-k1^2 + (h + k1) (h + k1\.b4 + k2)) (h + k2\.b4) + 
   k2 (-h k2\.b4 - k1 k2\.b4)), (
  h k2\.b4 + 
   k1 k2\.b4)/((-k1^2 + (h + k1) (h + k1\.b4 + k2)) (h + k2\.b4) + 
   k2 (-h k2\.b4 - k1 k2\.b4))}, {(
  k1 k2)/((-k1^2 + (h + k1) (h + k1\.b4 + k2)) (h + k2\.b4) + 
   k2 (-h k2\.b4 - k1 k2\.b4)), (
  h k2 + k1 k2)/((-k1^2 + (h + k1) (h + k1\.b4 + k2)) (h + k2\.b4) + 
   k2 (-h k2\.b4 - k1 k2\.b4)), (-k1^2 + (h + k1) (h + k1\.b4 + 
      k2))/((-k1^2 + (h + k1) (h + k1\.b4 + k2)) (h + k2\.b4) + 
   k2 (-h k2\.b4 - k1 k2\.b4))}}

So a solution with that fills a booklet written matrix element per matrix element.

So the best input for a lengthy output is now:

sol = DSolve[C2'[t] == h.C1[t] - B.C2'[t], C2[t], t]

So the best input for a teaching students instruction guide book is

sol = DSolve[D[C2m[t], t] == h.C1m[t] - Bm.D[C2m, t], C2m[t], t]


{{C2m[t] -> 
   C[1] + Inactive[Integrate][(-Bm.0 + h.C1m[K[1]]), {K[1], 1, t}]}}

Then substitute the matrix B and C1 into the solution.

One of these geniuses that already published a book about Mathematica figured the process, there is a in the real world to make the right reduces at the right time. To operate with symbolics is often better. To gain insight into tedious output is cumbersome.

In this question working with abbreviations can be fine in the pages long tedious solution. But that is still too long compared to not make use of the element per element matrix short again.

might be even better input style to get a more matrix-like solution.

MatrixExp[B*t]

{{RootSum[
   h^3 + h^2 k1 - h k1^2 + h^2 k1\.b4 + h k1 k1\.b4 + h^2 k2 + 
     h k1 k2 + h^2 k2\.b4 + h k1 k2\.b4 - k1^2 k2\.b4 + 
     h k1\.b4 k2\.b4 + k1 k1\.b4 k2\.b4 - 3 h^2 #1 - 2 h k1 #1 + 
     k1^2 #1 - 2 h k1\.b4 #1 - k1 k1\.b4 #1 - 2 h k2 #1 - k1 k2 #1 - 
     2 h k2\.b4 #1 - k1 k2\.b4 #1 - k1\.b4 k2\.b4 #1 + 3 h #1^2 + 
     k1 #1^2 + k1\.b4 #1^2 + k2 #1^2 + 
     k2\.b4 #1^2 - #1^3 &, (E^(t #1) h^2 + E^(t #1) h k1\.b4 + 
       E^(t #1) h k2 + E^(t #1) h k2\.b4 + E^(t #1) k1\.b4 k2\.b4 - 
       2 E^(t #1) h #1 - E^(t #1) k1\.b4 #1 - E^(t #1) k2 #1 - 
       E^(t #1) k2\.b4 #1 + E^(t #1) #1^2)/(3 h^2 + 2 h k1 - k1^2 + 
       2 h k1\.b4 + k1 k1\.b4 + 2 h k2 + k1 k2 + 2 h k2\.b4 + 
       k1 k2\.b4 + k1\.b4 k2\.b4 - 6 h #1 - 2 k1 #1 - 2 k1\.b4 #1 - 
       2 k2 #1 - 2 k2\.b4 #1 + 3 #1^2) &], -k1 RootSum[
    h^3 + h^2 k1 - h k1^2 + h^2 k1\.b4 + h k1 k1\.b4 + h^2 k2 + 
      h k1 k2 + h^2 k2\.b4 + h k1 k2\.b4 - k1^2 k2\.b4 + 
      h k1\.b4 k2\.b4 + k1 k1\.b4 k2\.b4 - 3 h^2 #1 - 2 h k1 #1 + 
      k1^2 #1 - 2 h k1\.b4 #1 - k1 k1\.b4 #1 - 2 h k2 #1 - k1 k2 #1 - 
      2 h k2\.b4 #1 - k1 k2\.b4 #1 - k1\.b4 k2\.b4 #1 + 3 h #1^2 + 
      k1 #1^2 + k1\.b4 #1^2 + k2 #1^2 + 
      k2\.b4 #1^2 - #1^3 &, (-E^(t #1) h - E^(t #1) k2\.b4 + 
        E^(t #1) #1)/(3 h^2 + 2 h k1 - k1^2 + 2 h k1\.b4 + 
        k1 k1\.b4 + 2 h k2 + k1 k2 + 2 h k2\.b4 + k1 k2\.b4 + 
        k1\.b4 k2\.b4 - 6 h #1 - 2 k1 #1 - 2 k1\.b4 #1 - 2 k2 #1 - 
        2 k2\.b4 #1 + 3 #1^2) &], 
  k1 k2\.b4 RootSum[
    h^3 + h^2 k1 - h k1^2 + h^2 k1\.b4 + h k1 k1\.b4 + h^2 k2 + 
      h k1 k2 + h^2 k2\.b4 + h k1 k2\.b4 - k1^2 k2\.b4 + 
      h k1\.b4 k2\.b4 + k1 k1\.b4 k2\.b4 - 3 h^2 #1 - 2 h k1 #1 + 
      k1^2 #1 - 2 h k1\.b4 #1 - k1 k1\.b4 #1 - 2 h k2 #1 - k1 k2 #1 - 
      2 h k2\.b4 #1 - k1 k2\.b4 #1 - k1\.b4 k2\.b4 #1 + 3 h #1^2 + 
      k1 #1^2 + k1\.b4 #1^2 + k2 #1^2 + k2\.b4 #1^2 - #1^3 &, 
    E^(t #1)/(3 h^2 + 2 h k1 - k1^2 + 2 h k1\.b4 + k1 k1\.b4 + 
        2 h k2 + k1 k2 + 2 h k2\.b4 + k1 k2\.b4 + k1\.b4 k2\.b4 - 
        6 h #1 - 2 k1 #1 - 2 k1\.b4 #1 - 2 k2 #1 - 2 k2\.b4 #1 + 
        3 #1^2) &]}, {-k1 RootSum[
    h^3 + h^2 k1 - h k1^2 + h^2 k1\.b4 + h k1 k1\.b4 + h^2 k2 + 
      h k1 k2 + h^2 k2\.b4 + h k1 k2\.b4 - k1^2 k2\.b4 + 
      h k1\.b4 k2\.b4 + k1 k1\.b4 k2\.b4 - 3 h^2 #1 - 2 h k1 #1 + 
      k1^2 #1 - 2 h k1\.b4 #1 - k1 k1\.b4 #1 - 2 h k2 #1 - k1 k2 #1 - 
      2 h k2\.b4 #1 - k1 k2\.b4 #1 - k1\.b4 k2\.b4 #1 + 3 h #1^2 + 
      k1 #1^2 + k1\.b4 #1^2 + k2 #1^2 + 
      k2\.b4 #1^2 - #1^3 &, (-E^(t #1) h - E^(t #1) k2\.b4 + 
        E^(t #1) #1)/(3 h^2 + 2 h k1 - k1^2 + 2 h k1\.b4 + 
        k1 k1\.b4 + 2 h k2 + k1 k2 + 2 h k2\.b4 + k1 k2\.b4 + 
        k1\.b4 k2\.b4 - 6 h #1 - 2 k1 #1 - 2 k1\.b4 #1 - 2 k2 #1 - 
        2 k2\.b4 #1 + 3 #1^2) &], 
  RootSum[h^3 + h^2 k1 - h k1^2 + h^2 k1\.b4 + h k1 k1\.b4 + h^2 k2 + 
     h k1 k2 + h^2 k2\.b4 + h k1 k2\.b4 - k1^2 k2\.b4 + 
     h k1\.b4 k2\.b4 + k1 k1\.b4 k2\.b4 - 3 h^2 #1 - 2 h k1 #1 + 
     k1^2 #1 - 2 h k1\.b4 #1 - k1 k1\.b4 #1 - 2 h k2 #1 - k1 k2 #1 - 
     2 h k2\.b4 #1 - k1 k2\.b4 #1 - k1\.b4 k2\.b4 #1 + 3 h #1^2 + 
     k1 #1^2 + k1\.b4 #1^2 + k2 #1^2 + 
     k2\.b4 #1^2 - #1^3 &, (E^(t #1) h^2 + E^(t #1) h k1 + 
       E^(t #1) h k2\.b4 + E^(t #1) k1 k2\.b4 - 2 E^(t #1) h #1 - 
       E^(t #1) k1 #1 - E^(t #1) k2\.b4 #1 + E^(t #1) #1^2)/(3 h^2 + 
       2 h k1 - k1^2 + 2 h k1\.b4 + k1 k1\.b4 + 2 h k2 + k1 k2 + 
       2 h k2\.b4 + k1 k2\.b4 + k1\.b4 k2\.b4 - 6 h #1 - 2 k1 #1 - 
       2 k1\.b4 #1 - 2 k2 #1 - 2 k2\.b4 #1 + 
       3 #1^2) &], -k2\.b4 RootSum[
    h^3 + h^2 k1 - h k1^2 + h^2 k1\.b4 + h k1 k1\.b4 + h^2 k2 + 
      h k1 k2 + h^2 k2\.b4 + h k1 k2\.b4 - k1^2 k2\.b4 + 
      h k1\.b4 k2\.b4 + k1 k1\.b4 k2\.b4 - 3 h^2 #1 - 2 h k1 #1 + 
      k1^2 #1 - 2 h k1\.b4 #1 - k1 k1\.b4 #1 - 2 h k2 #1 - k1 k2 #1 - 
      2 h k2\.b4 #1 - k1 k2\.b4 #1 - k1\.b4 k2\.b4 #1 + 3 h #1^2 + 
      k1 #1^2 + k1\.b4 #1^2 + k2 #1^2 + 
      k2\.b4 #1^2 - #1^3 &, (-E^(t #1) h - E^(t #1) k1 + 
        E^(t #1) #1)/(3 h^2 + 2 h k1 - k1^2 + 2 h k1\.b4 + 
        k1 k1\.b4 + 2 h k2 + k1 k2 + 2 h k2\.b4 + k1 k2\.b4 + 
        k1\.b4 k2\.b4 - 6 h #1 - 2 k1 #1 - 2 k1\.b4 #1 - 2 k2 #1 - 
        2 k2\.b4 #1 + 3 #1^2) &]}, {k1 k2 RootSum[
    h^3 + h^2 k1 - h k1^2 + h^2 k1\.b4 + h k1 k1\.b4 + h^2 k2 + 
      h k1 k2 + h^2 k2\.b4 + h k1 k2\.b4 - k1^2 k2\.b4 + 
      h k1\.b4 k2\.b4 + k1 k1\.b4 k2\.b4 - 3 h^2 #1 - 2 h k1 #1 + 
      k1^2 #1 - 2 h k1\.b4 #1 - k1 k1\.b4 #1 - 2 h k2 #1 - k1 k2 #1 - 
      2 h k2\.b4 #1 - k1 k2\.b4 #1 - k1\.b4 k2\.b4 #1 + 3 h #1^2 + 
      k1 #1^2 + k1\.b4 #1^2 + k2 #1^2 + k2\.b4 #1^2 - #1^3 &, 
    E^(t #1)/(3 h^2 + 2 h k1 - k1^2 + 2 h k1\.b4 + k1 k1\.b4 + 
        2 h k2 + k1 k2 + 2 h k2\.b4 + k1 k2\.b4 + k1\.b4 k2\.b4 - 
        6 h #1 - 2 k1 #1 - 2 k1\.b4 #1 - 2 k2 #1 - 2 k2\.b4 #1 + 
        3 #1^2) &], -k2 RootSum[
    h^3 + h^2 k1 - h k1^2 + h^2 k1\.b4 + h k1 k1\.b4 + h^2 k2 + 
      h k1 k2 + h^2 k2\.b4 + h k1 k2\.b4 - k1^2 k2\.b4 + 
      h k1\.b4 k2\.b4 + k1 k1\.b4 k2\.b4 - 3 h^2 #1 - 2 h k1 #1 + 
      k1^2 #1 - 2 h k1\.b4 #1 - k1 k1\.b4 #1 - 2 h k2 #1 - k1 k2 #1 - 
      2 h k2\.b4 #1 - k1 k2\.b4 #1 - k1\.b4 k2\.b4 #1 + 3 h #1^2 + 
      k1 #1^2 + k1\.b4 #1^2 + k2 #1^2 + 
      k2\.b4 #1^2 - #1^3 &, (-E^(t #1) h - E^(t #1) k1 + 
        E^(t #1) #1)/(3 h^2 + 2 h k1 - k1^2 + 2 h k1\.b4 + 
        k1 k1\.b4 + 2 h k2 + k1 k2 + 2 h k2\.b4 + k1 k2\.b4 + 
        k1\.b4 k2\.b4 - 6 h #1 - 2 k1 #1 - 2 k1\.b4 #1 - 2 k2 #1 - 
        2 k2\.b4 #1 + 3 #1^2) &], 
  RootSum[h^3 + h^2 k1 - h k1^2 + h^2 k1\.b4 + h k1 k1\.b4 + h^2 k2 + 
     h k1 k2 + h^2 k2\.b4 + h k1 k2\.b4 - k1^2 k2\.b4 + 
     h k1\.b4 k2\.b4 + k1 k1\.b4 k2\.b4 - 3 h^2 #1 - 2 h k1 #1 + 
     k1^2 #1 - 2 h k1\.b4 #1 - k1 k1\.b4 #1 - 2 h k2 #1 - k1 k2 #1 - 
     2 h k2\.b4 #1 - k1 k2\.b4 #1 - k1\.b4 k2\.b4 #1 + 3 h #1^2 + 
     k1 #1^2 + k1\.b4 #1^2 + k2 #1^2 + 
     k2\.b4 #1^2 - #1^3 &, (E^(t #1) h^2 + E^(t #1) h k1 - 
       E^(t #1) k1^2 + E^(t #1) h k1\.b4 + E^(t #1) k1 k1\.b4 + 
       E^(t #1) h k2 + E^(t #1) k1 k2 - 2 E^(t #1) h #1 - 
       E^(t #1) k1 #1 - E^(t #1) k1\.b4 #1 - E^(t #1) k2 #1 + 
       E^(t #1) #1^2)/(3 h^2 + 2 h k1 - k1^2 + 2 h k1\.b4 + 
       k1 k1\.b4 + 2 h k2 + k1 k2 + 2 h k2\.b4 + k1 k2\.b4 + 
       k1\.b4 k2\.b4 - 6 h #1 - 2 k1 #1 - 2 k1\.b4 #1 - 2 k2 #1 - 
       2 k2\.b4 #1 + 3 #1^2) &]}}

This already is longer than one page of letter size. IndentityMatrix[3]!

DSolve[C2'[t]==h.C1[t]-B.C2'[t],{CA2[t], CB2[t], CC2[t]},t]

{{CA2[t] -> 
   C[1] + Inactive[Integrate][
     h.h.{CB0 (-((E^(k1 K[1]) (h + k1) (h + k2))/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 
             3 h k1 k2)) + ((1 - E^((-h - k1) K[1])) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (h k2 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CC0 (((1 - E^((-h - k1) K[1])) k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k1 K[1]) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           h^2 + 2 h k1 + h k2 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CA0 (-((k1 k2)/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) - (
           E^(k1 K[1]) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 
            3 h k1 k2) + ((1 - E^((-h - k1) K[1])) (h^2 + h k1 + 
              2 h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)), 
       CB0 (((1 - E^((-h - k1 - k2) K[1])) (h + k1) (h + k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k1 K[1]) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k2 K[1]) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CC0 (-((E^(k1 K[1]) k1 k2)/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 
             3 h k1 k2)) + ((1 - E^((-h - k1 - k2) K[1])) (h k2 + 
              k1 k2))/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k2 K[1]) (h^2 + 2 h k1 + h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CA0 (-((E^(k2 K[1]) k1 k2)/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 
             3 h k1 k2)) + ((1 - E^((-h - k1 - k2) K[1])) (h k1 + 
              k1 k2))/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k1 K[1]) (h^2 + h k1 + 2 h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)), 
       CB0 (-((E^(k2 K[1]) (h + k1) (h + k2))/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) - (h k1 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 
            3 h k1 k2) + ((1 - E^((-h - k2) K[1])) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CC0 (-((k1 k2)/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) - (
           E^(k2 K[1]) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 
            3 h k1 k2) + ((1 - E^((-h - k2) K[1])) (h^2 + 2 h k1 + 
              h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CA0 (((1 - E^((-h - k2) K[1])) k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k2 K[1]) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           h^2 + h k1 + 2 h k2 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2))}/(
     1 + h), {K[1], 1, t}], 
  CB2[t] -> 
   C[2] + Inactive[Integrate][
     h.h.{CB0 (-((E^(k1 K[2]) (h + k1) (h + k2))/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 
             3 h k1 k2)) + ((1 - E^((-h - k1) K[2])) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (h k2 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CC0 (((1 - E^((-h - k1) K[2])) k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k1 K[2]) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           h^2 + 2 h k1 + h k2 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CA0 (-((k1 k2)/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) - (
           E^(k1 K[2]) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 
            
            3 h k1 k2) + ((1 - E^((-h - k1) K[2])) (h^2 + h k1 + 
              2 h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)), 
       CB0 (((1 - E^((-h - k1 - k2) K[2])) (h + k1) (h + k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k1 K[2]) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k2 K[2]) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CC0 (-((E^(k1 K[2]) k1 k2)/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 
             3 h k1 k2)) + ((1 - E^((-h - k1 - k2) K[2])) (h k2 + 
              k1 k2))/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k2 K[2]) (h^2 + 2 h k1 + h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CA0 (-((E^(k2 K[2]) k1 k2)/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 
             3 h k1 k2)) + ((1 - E^((-h - k1 - k2) K[2])) (h k1 + 
              k1 k2))/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k1 K[2]) (h^2 + h k1 + 2 h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)), 
       CB0 (-((E^(k2 K[2]) (h + k1) (h + k2))/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) - (h k1 + k1 k2)/(
           
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 
            3 h k1 k2) + ((1 - E^((-h - k2) K[2])) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CC0 (-((k1 k2)/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) - (
           E^(k2 K[2]) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 
            3 h k1 k2) + ((1 - E^((-h - k2) K[2])) (h^2 + 2 h k1 + 
              h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CA0 (((1 - E^((-h - k2) K[2])) k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k2 K[2]) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           h^2 + h k1 + 2 h k2 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2))}/(
     1 + h), {K[2], 1, t}], 
  CC2[t] -> 
   C[3] + Inactive[Integrate][
     h.h.{CB0 (-((E^(k1 K[3]) (h + k1) (h + k2))/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 
             3 h k1 k2)) + ((1 - E^((-h - k1) K[3])) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (h k2 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CC0 (((1 - E^((-h - k1) K[3])) k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k1 K[3]) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           h^2 + 2 h k1 + h k2 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CA0 (-((k1 k2)/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) - (
           E^(k1 K[3]) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 
            3 h k1 k2) + ((1 - E^((-h - k1) K[3])) (h^2 + h k1 + 
              2 h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)), 
       CB0 (((1 - E^((-h - k1 - k2) K[3])) (h + k1) (h + k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k1 K[3]) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k2 K[3]) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CC0 (-((E^(k1 K[3]) k1 k2)/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 
             3 h k1 k2)) + ((1 - E^((-h - k1 - k2) K[3])) (h k2 + 
              k1 k2))/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k2 K[3]) (h^2 + 2 h k1 + h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CA0 (-((E^(k2 K[3]) k1 k2)/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 
             3 h k1 k2)) + ((1 - E^((-h - k1 - k2) K[3])) (h k1 + 
              k1 k2))/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k1 K[3]) (h^2 + h k1 + 2 h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)), 
       CB0 (-((E^(k2 K[3]) (h + k1) (h + k2))/(
            h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) - (h k1 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 
            3 h k1 k2) + ((1 - E^((-h - k2) K[3])) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CC0 (-((k1 k2)/(h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) - (
           E^(k2 K[3]) (h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 
            3 h k1 k2) + ((1 - E^((-h - k2) K[3])) (h^2 + 2 h k1 + 
              h k2 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2)) + 
        CA0 (((1 - E^((-h - k2) K[3])) k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           E^(k2 K[3]) (h k1 + k1 k2))/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2) - (
           h^2 + h k1 + 2 h k2 + k1 k2)/(
           h^3 + 2 h^2 k1 + 2 h^2 k2 + 3 h k1 k2))}/(
     1 + h), {K[3], 1, t}]}}

Thanksss for all the help and insight. It seems to me that this is the answer im looking for. — felipe, Commented Jul 30, 2020 at 3:03

Stack Exchange Network

How to solve a Matrix of ODEs

1 Answer 1

Your Answer

Not the answer you're looking for? Browse other questions tagged
differential-equations
matrix
code-request
or ask your own question.

Hot Network Questions

How to solve a Matrix of ODEs

1 Answer 1

Your Answer

Sign up or log in

Post as a guest

Not the answer you're looking for? Browse other questions tagged differential-equationsmatrixcode-request or ask your own question.

Related

Hot Network Questions

Not the answer you're looking for? Browse other questions tagged
differential-equations
matrix
code-request
or ask your own question.