function iv = IV1()
  iv = zeros(11,1)
  iv( 0 +1) = 1
  iv(10 +1) = 4
endfunction

function iv = IV2()
  iv = zeros(11,1)
  iv( 0 +1) = 1
  iv(10 +1) = 4
  for x=1:9
    if x <= 5
      iv(x+1) = 2*x+1
    else
      iv(x+1) = -7/5 * x +18
    end
  end
endfunction

function u = ex1a(iv,time)
  
   dt = 0.25
   
   tmp = iv
   u = iv
   for t=1:time
    for i=1:9
      tmp(i + 1 ) = dt * ( u(i) -2 *u(i+1) + u(i+2)) + u(i+1)
    end
    u=tmp
  end
endfunction

function u = ex1b(iv,time)
   tmp = iv
   u = iv
   for t=1:time
    for i=1:9
      tmp(i + 1 ) = 0.5 * (u(i) + u(i+2))
    end
    u=tmp
  end
endfunction


function res = ex2(iv)
  n = size(iv,1)
  A = zeros(n,n)
  i=1
  A(1,1) = 1
  A(n,n) = 1
  for j=2:n-1
    A(j,i)   =  1
    A(j,i+1) = -2
    A(j,i+2) =  1
    i = i +1
  end

  printf("The arising matrix looks like this:\n")
  disp(A)

  // decomposing matrix into a,b,c
  a = zeros(n,1)  
  b = zeros(n,1)
  b(1) = A(1,1)
  b(n) = A(n,n)
  c = zeros(n,1)
  d = iv
  x = zeros(n,1)

  for i=2:n-1
    b(i) = A(i,i)
    a(i) = A(i,i-1)
    c(i) = A(i,i+1)
  end
  printf("a is\n")
  disp(a')
  printf("b is\n")
  disp(b')
  printf("c is\n")
  disp(c')
  
  // now compute new c, corresponding to the wikipedia entry
  for i=1:n-1
    if i == 1
      c(1) = c(1) / b(1)
    else
      c(i) = c(i) / (b(i) - c(i-1)*a(i))
    end
  end
  printf("new c is\n")
  disp(c')

  // now compute new d
  for i=1:n
    if i == 1
      d(1) = d(1) / b(1)
    else
      d(i) = (d(i) - d(i-1) *a(i)) / (b(i) - c(i-1)*a(i))
    end
  end
  printf("new d is\n")
  disp(d')

  // computing result
  for j=0:n-1
    i = n-j
    //printf("i=%d\n",i)
    if i == n
      x(n) = d(n)
    else
      x(i) = d(i) - c(i) *x(i+1)
    end
  end

  printf("\nown solution:\n")
  disp(x')
  printf("\ncorrect solution:\n")
  disp((A\iv)')
  res = x
endfunction




function plot1a()
  plot(ex1a(IV2(),0))                              
  plot(ex1a(IV2(),1))
  plot(ex1a(IV2(),10))
  plot(ex1a(IV2(),20))
endfunction


function plot1b()
  plot(ex1b(IV2(),0))                              
  plot(ex1b(IV2(),1))
  plot(ex1b(IV2(),10))
  plot(ex1b(IV2(),20))
endfunction

