I am tried to solve the following problem. here is the code and boundary conditions as well as parameters used in the problem. Please help me to get the numerical solution and getting plots between Cu and eta as well as D(f)(eta) vs eta.

restart;

Digits := trunc(evalhf(Digits));

15

ODEs := [diff(f(eta), `$`(eta, 3))+A^2+f(eta)*(diff(f(eta), `$`(eta, 2)))-(diff(f(eta), eta))^2+beta*((diff(g(eta), eta))^2-g(eta)*(diff(g(eta), `$`(eta, 2)))-1), lambda*(diff(g(eta), `$`(eta, 3)))+(diff(g(eta), `$`(eta, 2)))*f(eta)-g(eta)*(diff(f(eta), `$`(eta, 2)))];

`<,>`(ODEs[]);

Vector[column](%id = 18446744073898822582)

LB, UB := 0, 1;

BCs := [`~`[`=`](([D(f), f, g, (D@@2)(g)])(LB), [1+B1*((D@@2)(f))(0), 0, 0, 0])[], `~`[`=`](([D(f), D(g)])(UB), [A, 0])[]];

[D(f)(0) = 1 + B1 @@(D, 2)(f)(0), f(0) = 0, g(0) = 0,

@@(D, 2)(g)(0) = 0, D(f)(1) = A, D(g)(1) = 0]

Params := Record(A = .9, B1 = .5, beta = .5, lambda = .5);

NBVs := [-((D@@2)(f))(1) = `C*__f`];

Cf := `C*__f`;

Solve := module () local nbvs_rhs, Sol, ModuleApply, AccumData, ModuleLoad; export SavedData, Pos, Init; nbvs_rhs := `~`[rhs](:-NBVs); ModuleApply := subs(_Sys = {:-BCs[], :-NBVs[], :-ODEs[]}, proc ({ A::realcons := Params:-A, B1::realcons := Params:-B1, beta::realcons := Params:-beta, lambda::realcons := Params:-lambda }) Sol := dsolve(_Sys, _rest, numeric); AccumData(Sol, {_options}); Sol end proc); AccumData := proc (Sol::{Matrix, procedure, list({name, function} = procedure)}, params::(set(name = realcons))) local n, nbvs; if Sol::Matrix then nbvs := seq(n = Sol[2, 1][1, Pos(n)], n = nbvs_rhs) else nbvs := `~`[`=`](nbvs_rhs, eval(nbvs_rhs, Sol(:-LB)))[] end if; SavedData[params] := Record[packed](params[], nbvs) end proc; ModuleLoad := eval(Init); Init := proc () Pos := proc (n::name) local p; option remember; member(n, Sol[1, 1], 'p'); p end proc; SavedData := table(); return end proc; ModuleLoad() end module;

colseq := [red, green, blue, brown];

Pc := B1 = .5, A = .1, beta = .5;

Ps := [B1 = .5, A = .1, beta = .5];

Pv := [lambda = [.5, 1, 1.5, 2]];

for i to nops(Ps) do plots:-display([seq(plots:-odeplot(Solve(lhs(Pv[i]) = rhs(Pv[i])[j], Ps[i][], Pc), [eta, theta(eta)], 'color' = colseq[j], 'legend' = [lhs(Pv[i]) = rhs(Pv[i])[j]]), j = 1 .. nops(rhs(Pv[i])))], 'axes' = 'boxed', 'gridlines' = false, 'labelfont' = ['TIMES', 'BOLDOBLIQUE', 16], 'caption' = nprintf(cat(`$`("\n%a = %4.2f, ", nops(Ps[i])-1), "%a = %4.2f\n\n"), `~`[lhs, rhs](Ps[i])[]), 'captionfont' = ['TIMES', 16]) end do;

Error, (in dsolve/numeric/process_input) invalid argument: (B1 = .5)[]

Please help me to get the graph of CU v/s eta also D(f)(eta) vs eta