title: Solving systems of second order differential equations
body: <p>I'm working on a script in mathematica that will take simulate a string held at either end and plucked, by solving the wave equation via numerical methods. (<a href="http://en.wikipedia.org/wiki/Wave_equation#Investigation_by_numerical_methods" rel="nofollow">http://en.wikipedia.org/wiki/Wave_equation#Investigation_by_numerical_methods</a>) </p>
  
  <pre><code>n = 5; (*The number of discreet elements to be used*)
  L = 1.0;  (*The length of the string that is vibrating*)
  a = 1.0/3.0; (*The distance from the left side that the string is \
  plucked at*)
   T = 1; (*The tension in the string*)
  [Rho] = 1; (*The length density of the string*)
  y0 = 0.1; (*The vertical distance of the string pluck*)
  [CapitalDelta]x = L/n; (*The length of each discreet element*)
  m = ([Rho]*L)/n;(*The mass of each individual node*)
  c = Sqrt[T/[Rho]];(*The speed at which waves in the string propogate*)
  </code></pre>
  
  <p>I set all my variables</p>
  
  <pre><code>Y[t] = Array[f[t], {n - 1, 1}];
  MatrixForm(*Creates a vector size n-1 by 1 of functions \
  representing each node*)
  </code></pre>
  
  <p>I define my Vector of nodal position functions</p>
  
  <pre><code>K = MatrixForm[
  SparseArray[{Band[{1, 1}] -&gt; -2, Band[{2, 1}] -&gt; 1, 
  Band[{1, 2}] -&gt; 1}, {n - 1, 
  n - 1}]](*Creates a matrix size n by n governing the coupling \
  between each node*)
  </code></pre>
  
  <p>I create the stiffness matrix relating all the nodal functions to one another</p>
  
  <pre><code>Y0 = MatrixForm[
  Table[Piecewise[{{(((i*L)/n)*y0)/a, 
    0 &lt; ((i*L)/n) &lt; a}, {(-((i*L)/n)*y0)/(L - a) + (y0*L)/(L - a), 
    a &lt; ((i*L)/n) &lt; L}}], {i, 1, n - 1}]]
  </code></pre>
  
  <p>I define the initial positions of each node using a piecewise function</p>
  
  <pre><code>NDSolve[{Y''[t] == (c/[CapitalDelta]x)^2 Y[t].K, Y[0] == Y0, 
  Y'[0] == 0}, 
  Y, {t, 0, 10}];(*Numerically solves the system of second order DE's*)
  </code></pre>
  
  <p>Finally, This should solve for the values of the individual nodes, but it returns an error:</p>
  
  <p>"NDSolve::ndinnt : Initial condition [Y0 table] is not a number or a rectangular array"</p>
  
  <p>So , it would seem that I don't have a firm grasp on how matrices work in mathematica.  I would greatly appreciate it if anyone could help me get this last line of code to run properly.</p>
  
  <p>Thank you,
  Brad</p>