restart;
#Setup procedure for Midpoint method for solving ode: dy/dx=y
#The procedure accepts as inputs: n, the number of steps, h, the stepsize, and x0,y0, the initial condition
MidPoint := proc (n,h,x0,y0)
x[1] := x0; y[1] := y0;
printf("[%4.4f",x[1]); #these "printf" statements just give me a way of outputting the data as I like it
for i from 2 to n do
x[i]:=x[i-1]+h; #here we set up all the x values
printf( ",%4.4f",x[i]);
end do;
printf("]\134n");
printf("[%4.4f",y[1]);
for j from 2 to n do
k1:=(h/2)*y[j-1]; #here we have "y[i-1]" because that r.h.s. of our ode is "y"
y[j] :=y[j-1] + h*(y[j-1]+k1); #again, f(x[i-1]+(h/2),y[i-1]+k)=y[i-1]+k
printf(",%4.4f", y[j]);
end do;
printf("]\134n");
end proc:
Warning, `x` is implicitly declared local to procedure `MidPoint`
Warning, `y` is implicitly declared local to procedure `MidPoint`
Warning, `i` is implicitly declared local to procedure `MidPoint`
Warning, `j` is implicitly declared local to procedure `MidPoint`
Warning, `k1` is implicitly declared local to procedure `MidPoint`
MidPoint(4,0.1,0,1);
[0.0000,0.1000,0.2000,0.3000]
[1.0000,1.1050,1.2210,1.3492]
with(plots):
Xvalues:=[0.0000,0.1000,0.2000,0.3000];Yvalues:=[1.0000,1.1050,1.2210,1.3492];
NyYkIiIhRiQkIiUrNSEiJSQiJSs/RickIiUrSUYn
NyYkIiYrKyIhIiUkIiZdNSJGJSQiJjVBIkYlJCImI1w4RiU=
points:= {seq([Xvalues[n],Yvalues[n]],n=1..4)}:
numeric_sol:=pointplot(points):
#The exact solution is
with(DEtools):
DE:=diff(y(x),x)-y(x);
IC:=y(0)=1;
LCYtSSVkaWZmRyUqcHJvdGVjdGVkRzYkLUkieUc2IjYjSSJ4R0YpRisiIiJGJyEiIg==
Ly1JInlHNiI2IyIiISIiIg==
dsolve({DE,IC});
Ly1JInlHNiI2I0kieEdGJS1JJGV4cEc2JCUqcHJvdGVjdGVkR0koX3N5c2xpYkdGJUYm
exact_sol:=plot(exp(x),x=0..0.3):
display({exact_sol,numeric_sol})
Warning, inserted missing semicolon at end of statement
NiktJSdQT0lOVFNHNiY3JCQiIiEhIiIkIiM1ISIiNyQkIiIiISIiJCIlMDYhIiQ3JCQiIiMhIiIkIiVANyEiJDckJCIiJCEiIiQiJiNcOCEiJS0lJ0NVUlZFU0c2JDdTNyQkIiIhISIiJCIjNSEiIjckJCIxKioqKioqXGk5UmwhIz0kIjE4Zjl0MGMxNSEjOjckJCItdlZBKUdBIiEjOCQiMiNvP0ArUkk3NSEjOzckJCIsdiRldWk9ISM3JCIyYnpYPC4tKT01ISM7NyQkIixEJzMmb10jISM3JCIyZnFmZU8mUUQ1ISM7NyQkIjInKioqXChvWCp5OSQhIz0kIjE8XHVaJ3o+LiIhIzo3JCQiMicqKipcUDlDQXUkISM9JCIyKCoqek1xNzhRNSEjOzckJCIyJSoqKlxQKnpoZFYhIz0kIjI7KilSVWNSWC8iISM7NyQkIjIlKioqXFA+ZlMqXCEjPSQiMlU3X1VrMzcwIiEjOzckJCIyJyoqKlwoPSRmJUdjISM9JCIyZHRsO3IpKnkwIiEjOzckJCIrRHksIkcnISM2JCIyS05mSnBDWzEiISM7NyQkIjEqKioqXDc8emJvISM8JCIxJFtLLWxpNDIiISM6NyQkIit2NCZHXSghIzYkIjI7T0k6KVsiejIiISM7NyQkIjEqKioqKlw3bkQ6KSEjPCQiMmQjKj4hcDUlXDMiISM7NyQkIitEISpveSgpISM2JCIxJ1syTlZiPDQiISM6NyQkIjEpKipcUHBuc00qISM8JCIyWyIpUUZnISl6NCIhIzs3JCQiLERGT0IrIiEjNiQiMmZUIkdZIkhhNSIhIzs3JCQiMioqKioqKioqUTUnZjUhIzwkIjIoUiNlJ2YmeTw2IiEjOzckJCIyKioqXFAvUUJFNiEjPCQiMjtYYnFMNSM+NiEjOzckJCIyKCoqKioqXCJvPyY9IiEjPCQiMjhpPHg6SWU3IiEjOzckJCIudlZiNCpcNyEjOCQiMidvNkcvI1FKOCIhIzs3JCQiLkRKJz1fNjghIzgkIjJ4KG8rI0dULDkiISM7NyQkIi12VnkhZVAiISM3JCIyPCFmJT1TJVxaNiEjOzckJCIudj1XVVtWIiEjOCQiMlVxM1JpKUdhNiEjOzckJCIyKSoqKlw3Qj4mKVwiISM8JCIxYThvOUFtaDYhIzo3JCQiLnYkPjptazohIzgkIjJpTisyPXIkcDYhIzs3JCQiLkRjZFFBaSIhIzgkIjJ4OCQqZlxCaDwiISM7NyQkIi12dExVJW8iISM3JCIyQjVvTSkqZk09IiEjOzckJCIyKSoqKioqXE5tJ1s8ISM8JCIyempcc04oMyI+IiEjOzckJCIyKSoqKlwoeWJeNj0hIzwkIjInM2JKQW9mKT4iISM7NyQkIjIpKipcUE1hS3M9ISM8JCIxJG9qJHl3IWY/IiEjOjckJCItRDZXJSlSPiEjNyQiMjgpMzgkUnhTQCIhIzs3JCQiMicqKioqKlxAODArIyEjPCQiMlZQISpSV2w5QSIhIzs3JCQiMikqKioqXDcsSGw/ISM8JCIxKDRILVMuJUg3ISM6NyQkIjInKipcUDR3KVI3IyEjPCQiMllBaDMiNGtPNyEjOzckJCIsdlpmIik9IyEjNiQiMmpPdFA9LVlDIiEjOzckJCIyKCoqXFAvLWFbQSEjPCQiMmJeN2ciKlJARCIhIzs3JCQiLnY9WWI7SiMhIzgkIjFCd2dAeTFnNyEjOjckJCIyJioqKipcaUBPdEIhIzwkIjJfbT5QS255RSIhIzs3JCQiMigqKlxQZkwnelYjISM8JCIybiRwWzZXM3c3ISM7NyQkIioqPj0rRCEiKiQiMVEqeVt5W1NHIiEjOjckJCItREUmNFFjIyEjNyQiMjsicEY9XEMjSCIhIzs3JCQiLnZWPjVwaSMhIzgkIjJ5SG8zJltVKzghIzs3JCQiK2JKKltvIyEjNSQiMml6bzwrKCl6SSIhIzs3JCQiLURyIls4diMhIzckIjJXNzh4PDNuSiIhIzs3JCQiMigqKioqKipIank1RyEjPCQiMiZcLXpgeGJDOCEjOzckJCIudlYhKWZUKEchIzgkIjIqeSs9ImZ5SEwiISM7NyQkIi5EY0k7WyRIISM4JCIyWXJUVGEpM1Q4ISM7NyQkIiIkISIiJCIySyt3diEpZSlcOCEjOy0lJkNPTE9SRzYmJSRSR0JHJCIjNSEiIiQiIiEhIiIkIiIhISIiLSUlVklFV0c2JDskIiIhISIiJCIiJCEiIiUoREVGQVVMVEctJStBWEVTTEFCRUxTRzYnSSJ4RzYiUSE2Ii0lJUZPTlRHNiUlIUclIUciIzUlK0hPUklaT05UQUxHJStIT1JJWk9OVEFMRy0lKkFYRVNTVFlMRUc2IyUnTk9STUFMRy0lKFNDQUxJTkdHNiMlLlVOQ09OU1RSQUlORURHLSUlUk9PVEc2Jy0lKUJPVU5EU19YRzYjJCIkXSQhIiItJSlCT1VORFNfWUc2IyQiIyEqISIiLSUtQk9VTkRTX1dJRFRIRzYjJCIlK04hIiItJS5CT1VORFNfSEVJR0hURzYjJCIlXU0hIiItJSlDSElMRFJFTkc2Ig==