(* 
      Introduction to using Mathematica for 4  
        ME6 Stress Analysis *) 
 (* 
      In the derivations I use a very small part of what Mathematica can do. *) 
 (* 
      You will have guessed that parenthesis - 
        asterisk is used for comments *)  
 (* Mathematica understands about most functions .  
To run a CELL,  
      you need to press Shitf + Enter  -  
        top status bar should flash 
' RUNNING',  
      and come up with a plot below *) 
Plot[{Sin[x], Cos[x], Exp[x], Cosh[x]}, {x, 0, 
        Pi}]

Plot[{Sin[x], Cos[x], 
        Exp[x], Cosh[x]}, {x, 0, Pi}, 
      PlotStyle -> {Hue[0.1], Hue[0.3], Hue[0.5], 
          Hue[0.7]}]

(* 
      Mathematica can also build lists using Table: *) 
Table[Sin[n Pi/6], {n, 0, 6}]
Plot[{Sin[x], Sin[2  x], Sin[3  x], Sin[4  x]}, {x, 0, Pi}, 
      PlotStyle -> Table[Hue[0.25 j], {j, 1, 4}]]

(* 
      It can also build arrays: *) 
consts = Array[c, 5]

u = 
        Sin[rho] (3 + Cos[th]);
v = Cos[rho] (3 + Sin[th]);
z = Sin[th];
ParametricPlot3D[{u, v, z}, {rho, 0, 2  Pi}, {th, 0, 
        2  Pi}]

D[Sin[x], x]
D[Sin[x], {x, 2}]
Integrate[Sin[x], x]
Integrate[Sin[x], {x, 0, Pi}]

Series[
      Tan[x], {x, 0, 15}]

Solve[x^2 + 5  x - 6 == 0, x]

(* 
      It can also solve differential equations! *) 
DSolve[u''[x] - 9  u[x] == 0, u[x], 
      x]

DSolve[u''[x] - 4  u'[x] + 4  u[x] == 0, u[x], 
      x]

(* Mathematica is very good at solving systems of equations. Command SOLVE needs a list of equations, followed by a list of unknowns to find. *)
SOLUTION = 
        Solve[{x + 3  y + z == 2, 
-x - 2  y + z == 
              5, 
3  x + 7  y + z == -3}, {x, y, 
            z}]

(* 
      We call the result ' 
        SOLUTION' .  It contains RULES .  
          Each RULE describes a substitution :  
          it shows,  using an arrow ->,  
      which values should be given to each unknown . 
Our SOLUTION contains the rules enclosed in double curly brackets .  To extract the rules,  
      we needto use double square brackets,  
      and specify we need the first element *)  
SOLUTION[[1]]

{x, y, z} /. 
      SOLUTION[[1]]

OddSum[
        n_] := 
Module[{i, 
          total = 0}, 
Do[
          If[OddQ[i], total += i, 0], 
{i, 
            n}]; 
total
]

OddSum[100]