How to solve differential equation where constant is a changing index

Hello,

I have never used an ODE solver before and I'm having a fair bit of trouble. I have a differential equation that looks like:

, where , ρ, and g are all defined constants.

The issue I think I'm running into is that b is a vector, a function of x. At every point that differential equation is evaluated, the value of s and need to be the value at that particular index. More generally, the differential equations looks something like this:

and needs to be evaluted iteratively for each element/index within , where the output is a vector the size of . At least that how I believe it should be solved.

I have tried using for loops and following the "ODE with time dependent terms" but I haven't had any success. Here is my last attempt at code if it helps.

Thanks!

function dsdx = myode(t,s,bed) %this function is a separate file
bed = interp1(x,bed,t)
dsdx = (tau_y / p*g.*(s-bed))
[t,s] = ode45(@(t,s) myode(t,s,x,bed),[0 100],s0);
end
======================================
clear
close all

a = 1000; %length of bed
b = 200; %height of bed
x = linspace(0,1000,100);

bed = ((-.2*x + b) + (a/100.*rand(1,1).*sin(x)) + (a/100.*rand(1,1).*cos(x)))./(exp(x./a)); %Generate random bed terrain
bed = bed + (2.*rand(1,numel(bed))-1).*bed./10; %Add noise to generate  random bed terrain
bed = bed + (2.*rand(1,numel(bed))-1).*bed./100; %Add further noise

x = x';
bed = bed';

tau_y = 150*1000; %kPa
p = .917; %kg/m^3
g = 9.81; %m/s^2

s0 = bed(1);

[t,s] = ode45(@(t,s) myode(t,s,x,bed),x,s0);

[Walter Roberson]

function dsdx = myode(t,s,bed) %this function is a separate file

Okay, you have myode.m

[t,s] = ode45(@(t,s) myode(t,s,x,bed),[0 100],s0);

and inside it, you call ode45, asking ode45 to process the function myode... which is the same function you are already in.

Answer 1

[Torsten]

rng('default')
a = 1000; %length of bed
b = 200; %height of bed
x = linspace(0,1000,100); 
bed = ((-.2*x + b) + (a/100.*rand(1,1).*sin(x)) + (a/100.*rand(1,1).*cos(x)))./(exp(x./a)); %Generate random bed terrain
bed = bed + (2.*rand(1,numel(bed))-1).*bed./10; %Add noise to generate  random bed terrain
bed = bed + (2.*rand(1,numel(bed))-1).*bed./100; %Add further noise
x = x';
bed = bed';
tau_y = 150*1000; %kPa
p = .917; %kg/m^3
g = 9.81; %m/s^2
s0 = bed(1);
[t,s] = ode15s(@(t,s) myode(t,s,x,bed,tau_y,p,g),x,s0);
plot(t,s)

function dsdx = myode(t,s,x,bed,tau_y,p,g) %this function is a separate file
  bedt = interp1(x,bed,t);
  dsdx = tau_y / (p*g*(s-bedt));
end

[Tanner May]

This definitely helped steer me in the right direction. Thanks!