function [] = convolveToSolve2250()
clc % this clears the command line

n = 10000;
%define time row vector
t = linspace(0,30,n);

%% define forcing vector
b1 = 1; % turn-on time
b2 = 1; % amplitude
k = 8; %steepness of logistic function
f = forcing(t,b1,b2,k);


%% define kernel
a = -1/5;
omega = 3*pi;
k = kernel(t,a,omega);


%% convolve the kernel with the forcing:
x = convolver(t,f,k);
size(t)
size(x)


%%% plot everything
figure(1)
subplot(3,1,1)
plot(t,k,'k','linewidth',2)
set(gca,'fontsize',16)
xlabel('\tau')
ylabel('Kernel of DE')

subplot(3,1,2)
plot(t,f,'k','linewidth',2)
set(gca,'fontsize',16)
xlabel('t')
ylabel('Forcing f(t)')
ylim([-.2 1.2])

subplot(3,1,3)
plot(t,x,'k','linewidth',2)
set(gca,'fontsize',16)
xlabel('t')
ylabel('DE solution x(t)')
ylim([-.2 .3])









%%%%%%%%%%%% define sub-functions

%% forcing:
function [f]  = forcing(t,b1,b2,k)
    f = b2.*heaviside(t-b1+.000000003); %step function
    %f=b2.*(1./(1+exp(-2*k.*(t-b1)))); %logistic function
    f = double(f);
    
%% kernel:
function [k] = kernel(t,a,omega)
    k = exp(a*t).*sin(omega*t);
    
function [f_star_k] = convolver(t,f,k)
    dt = t(2)-t(1);
    f_star_k = conv(f,k)*dt;
    f_star_k = f_star_k(1:length(t));
     

    