matlab:打靶法求解边值问题

函数1：龙格库塔方法
%书籍：常用数值算法及其matlab实现
%第10章 常微分方程边值问题的数值解法,例11.1使用
%打靶法+四阶龙格库塔方法
%调用：function [t,z] = rk4symeq(fun, t0, tf, Za, h)
function [t,w] = shootin(fun1, fun2,t0, tf, alpha, beta, h)
%fun1, fun2:边值问题转化为初值问题的函数
%t0, tf：区间
%alpha, beta,初值和终值,是一个数
%h：步长
M = floor(tf-t0)/h ;      %离散点的个数M+1
if t0 >= tf
    printf('左端点必须小于右端点');
    return;
end
N = 2;           %获得变量个数,N =2 ,二阶方程
z1 = zeros(M+1, N);
z2 = zeros(M+1, N);
t = zeros(M +1, 1);
t =[t0 : h :tf]';
z1(1,:) = [alpha; 0] ;            %假设alpha为数字，补充导数的初值，变成列向量，符合rk4symeq的要求
z2(1,:) = [0; 1 ];            %假设alpha为数字，补充导数的初值，变成列向量，符合rk4symeq的要求
[t,z1] = rk4symeq(fun1, t0, tf, [alpha; 0], h);  %第1个微分方程
[t,z2] = rk4symeq(fun2, t0, tf, [0; 1], h);   %第2个微分方程

U = z1(:,1);        %取出y1，是一个列向量
V = z2(:,1);        %取出y2，是一个列向量
w = U + (