2
$\begingroup$

I would like to construct a matrix $G$ composed of block matrices $G^{(i)}$ as defined entry-wise below stacked vertically as $$G = \begin{bmatrix} G^{(1)} \\ G^{(2)} \\ \vdots \\ G^{(i)} \\ \vdots \\ G^{(n)} \end{bmatrix}$$ for some natural number $n$.

$G^{(i)}_{j,k}$ stands for the element of the $j$'th row and the $k$'th column of the $i$'th block matrix. \begin{align} G^{(i)}_{j,k} :=& \frac{t_{i+1}-t}{h_i}\delta_{k,i}+\frac{t-t_i}{h_i}\delta_{k,i+1} -\frac16(t-t_i)(t_{i+1}-t)\Big(1+\frac{t_{i+1}-t}{h_i}\Big)\delta_{k,n+i} \\ &-\frac16(t-t_i)(t_{i+1}-t)\Big(1+\frac{t-t_i}{h_i}\Big)\delta_{k,n+i+1} \ \Bigg|_{t=t_{j,i}} \end{align} where $i\in\{1,2,\cdots,n\}$, $j\in\{1,2,\cdots,m(j)\}$, $m(j)$ is a natural numbered function of $j$, and $k\in\{1,2,\cdots,2n-2\}$, $\delta_{i,j}$ is the Kronecker delta function, and the vertical bar at the right end indicates $t$ is to be evaluated at $t_{j,t}$.

What is the most elegant and convenient way to do this?

$\endgroup$
6
  • $\begingroup$ ` t is to be evaluated at t_j,t` - explain please. You have a mixture of $t,t_i,t_{i+1}$ in your problem, but is $t$ a function? And if so why does it take either 1 (i or i+1) or two arguments (j,i)? $\endgroup$
    – flinty
    Commented Aug 4, 2020 at 17:49
  • $\begingroup$ @flinty: It just means substituting $t$ by $t_{j,i}$. $\endgroup$
    – Hans
    Commented Aug 4, 2020 at 19:23
  • $\begingroup$ Use SparseArray. Or DiagonalMatrix. What are the limits on j and k? Are the G's square matrices? Can n be anything or is it related to the bounds on i and j? More details would be helpful. $\endgroup$
    – march
    Commented Aug 4, 2020 at 21:08
  • $\begingroup$ You have a $\delta_{k,i}$ and a $\delta_{n+k,n+i}$, but $\delta_{n+k,n+i}=\delta_{k,i}$, so these actually are part of the same matrix element. I suspect that you meant something else, however, since you wrote them separately. Can you clarify? I think editing your post with a small example output would be very helpful, with explicit (small) limits on j and j like 3 or something, i = 2 (so, two blocks), and n = something. As it stands, $n$ actually doesn't do anything in this expression, as far as I can tell. $\endgroup$
    – march
    Commented Aug 4, 2020 at 21:19
  • $\begingroup$ @march: You are absolutely right. I have edited my question. Please review. Thank you. $\endgroup$
    – Hans
    Commented Aug 4, 2020 at 23:18

2 Answers 2

3
$\begingroup$

Assuming you have definitions already for n, single argument t[x], and two-argument t[x,y], then this is fairly straightforward:

g[i_, n_] := 
 Table[(t[i + 1] - t[j, i])/h[i] KroneckerDelta[k, 
     i] + (t[j, i] - t[i])/h[i] KroneckerDelta[k, i + 1] - 
   1/6 (t[j, i] - t[i]) (t[i + 1] - 
      t[j, i]) (1 + (t[i + 1] - t[j, i])/h[i]) KroneckerDelta[k, 
     n + i] - 
   1/6 (t[j, i] - t[i]) (t[i + 1] - 
      t[j, i]) (1 + ((t[j, i] - t[i])/h[i])) KroneckerDelta[k, 
     n + i + 1], {j, 1, 3}, {k, 1, 2*n-2}]

result = With[{n = 4}, Join @@ Array[g[#, n] &, n]]
result // Dimensions (*expected {12,3}*)
result // MatrixForm
$\endgroup$
12
  • $\begingroup$ Is there a command to show the final explicit matrix form without the explicit KroneckerDelta function? $\endgroup$
    – Hans
    Commented Aug 4, 2020 at 19:32
  • 1
    $\begingroup$ @Hans can you explain what you're after? MatrixForm is just a presentational thing - you shouldn't use it in computations. $\endgroup$
    – flinty
    Commented Aug 4, 2020 at 22:05
  • $\begingroup$ Your solution is good. I understand I can use it for computation. I am wondering if there is a command to see the matrix in the explicit more conventional form as a sanity check. $\endgroup$
    – Hans
    Commented Aug 4, 2020 at 22:55
  • $\begingroup$ Yes that's what the //MatrixForm is for. It just looks weird because the expressions in the matrix are very wide. $\endgroup$
    – flinty
    Commented Aug 4, 2020 at 22:57
  • 1
    $\begingroup$ @Hans Array[g[#, n] &, n] generates a list of {g[1,n],g[2,n],g[3,n],..., g[n,n]} and @@ splices this list into a sequence of arguments for Join , so Join[g[1,n],g[2,n],g[3,n]...]. This is stacking them 'vertically' though Mathematica isn't like matlab; it works on lists, so it doesn't have a notion of vertical/horizontal, just lists and list nesting. If you wanted them horizontal you would use Join[g[1,n],g[2,n],g[3,n]....,{2}] with a level spec {2} at the end. $\endgroup$
    – flinty
    Commented Aug 5, 2020 at 11:26
2
$\begingroup$

Is this a formula with the Einstein sum convention. If so then the summation over two same indices has to be carried through.

g[i_, j_, k_, n_] := 
 Sum[(t[i + 1] - t[j, i])/h[ii] KroneckerDelta[k, 
     ii] + (t[j, i] - t[ii])/h[ii] KroneckerDelta[k, ii + 1] - 
   1/6 (t[j, i] - t[ii]) (t[ii + 1] - 
      t[j, i]) (1 + (t[ii + 1] - t[j, i])/h[ii]) KroneckerDelta[n + k,
      n + ii] - 
   1/6 (t[j, i] - t[ii]) (t[ii + 1] - 
      t[j, i]) (1 + ((t[j, i] - t[ii])/h[ii])) KroneckerDelta[n + k, 
     n + ii + 1], {ii, 1, 3}]

g[1, j, k, n]

(KroneckerDelta[1, k] (t[2] - t[j, 1]))/h[1] + (
 KroneckerDelta[2, k] (t[2] - t[j, 1]))/h[2] + (
 KroneckerDelta[3, k] (t[2] - t[j, 1]))/h[3] + (
 KroneckerDelta[2, k] (-t[1] + t[j, 1]))/h[1] - 
 1/6 KroneckerDelta[1 + n, 
   k + n] (1 + (t[2] - t[j, 1])/h[1]) (t[2] - t[j, 1]) (-t[1] + 
    t[j, 1]) + (KroneckerDelta[3, k] (-t[2] + t[j, 1]))/h[2] - 
 1/6 KroneckerDelta[2 + n, 
   k + n] (1 + (t[3] - t[j, 1])/h[2]) (t[3] - t[j, 1]) (-t[2] + 
    t[j, 1]) + (KroneckerDelta[4, k] (-t[3] + t[j, 1]))/h[3] - 
 1/6 KroneckerDelta[3 + n, 
   k + n] (1 + (t[4] - t[j, 1])/h[3]) (t[4] - t[j, 1]) (-t[3] + 
    t[j, 1]) - 
 1/6 KroneckerDelta[2 + n, 
   k + n] (t[2] - t[j, 1]) (-t[1] + t[j, 1]) (1 + (-t[1] + t[j, 1])/
    h[1]) - 1/
  6 KroneckerDelta[3 + n, 
   k + n] (t[3] - t[j, 1]) (-t[2] + t[j, 1]) (1 + (-t[2] + t[j, 1])/
    h[2]) - 1/
  6 KroneckerDelta[4 + n, 
   k + n] (t[4] - t[j, 1]) (-t[3] + t[j, 1]) (1 + (-t[3] + t[j, 1])/
    h[3])

The interpretation distinguishes between the indices i of the Einstein summation and the indices i from the discretization of time i in t_ij.

Example for displaying this as a matrix is

Table[h[i, j, k, n], {i, 1, 3}] // MatrixForm

enter image description here

This is much to lengthy for g.

$\endgroup$
1
  • $\begingroup$ Actually, there is no summation over 'ii' here. 'ii' should be the index of the block matrix. $\endgroup$
    – Hans
    Commented Aug 4, 2020 at 19:29

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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