I need help to plot this equation with Bvp4c

I just wanted to plot this entropy generation number with bvp4c. But I don't know how. Kindly help me

My code is

K = 0.1; We = 0.5; beta = 0.1; n=1; M=0.5;lambda=0.1;
B=0.5;Rd=0.7;Pr=3;Nb=0.5;Nt=0.7;Ec=0.8;Sc=1;m=0.5;B_1=10;B_2=10;
Gamma=0.1; theta_f=1.5;A=0.5;Br = 0.5; Re = 1; omega_e = 0.7; omega_c = 0.5;
infinity = 15;
maxinfinity = 20;
solinit = bvpinit(linspace(0,infinity,300),[1 1 0 0 1 0 0 0 1 0 1 0]);
sol = bvpsolver(@fsode,@fsbc,solinit);
eta = sol.x;
y = sol.y;

plot(eta,y(9,:),'g',eta(end),y(9,end));
legend('B = 0.5','B = 1.0','B = 1.5','B = 2.0',14);     
hold on

for Bnew = infinity+1:maxinfinity
    solinit = bvpxtend(sol,Bnew); 
    sol = bvp4c(@nanoode,@nanobc,solinit);
eta = sol.x; 
y = sol.y;
 end
 
 function dfdeta = fsode(eta,y)
    dfdeta = [  y(2)
                y(3)
              (((1+K)*(1+(1/beta))*(1+(1-n)*(We*y(3)^n)))^(-1))*(-(y(1)+y(4))*((1+(We*y(3))^n)^2)+(y(2)*y(2)+lambda*(y(2)+(eta/2)*y(3))+K*y(8)+M*y(2))*((1+(We*y(3))^n)^2))
                y(5)
                y(6)
              (((1+K)*(1+(1/beta))*(1+(1-n)*(We*y(6)^n)))^(-1))*(-(y(1)+y(4))*((1+(We*y(6))^n)^2)+(y(5)*y(5)+lambda*(y(5)+(eta/2)*y(6))+M*y(5))*((1+(We*y(6))^n)^2))
                y(8)
                ((1+K/2)^(-1))*(-(y(1)+y(4))*y(8)+y(7)*y(2)+(1/2)*(3*y(7)+eta*y(8))+K*B*(2*y(7)-y(3)))
                y(10)
                (Pr*((y(1)+y(4))*y(10)-(eta/2)*A*y(10)+Nb*y(10)*y(12)+Nt*y(10)*y(10))-Pr*(2+K)*(1+1/beta)*Ec*(y(3)*y(3)+y(6)*y(6))-Pr*M*Ec*(y(2)*y(2)+y(5)*y(5))-(1/2)*Pr*K*Ec*(y(3)-2*y(7))*(y(3)-2*y(7)))/(1+Rd*((theta_f-1)*y(9))^3)
                y(12)
                -(Nt/Nb)*((Pr*((y(1)+y(4))*y(10)-(eta/2)*A*y(10)+Nb*y(10)*y(12)+Nt*y(10)*y(10))-Pr*(2+K)*(1+1/beta)*Ec*(y(3)*y(3)+y(6)*y(6))-Pr*M*Ec*(y(2)*y(2)+y(5)*y(5))-(1/2)*Pr*K*Ec*(y(3)-2*y(7))*(y(3)-2*y(7)))/(1+Rd*((theta_f-1)*y(9))^3))+Sc*Pr*((y(1)+y(4))-(eta/2)*A)*y(12)+K*y(11)
                ];
                
  end
 function res = fsbc(y0,yinf)
    res = [ y0(1)-1
            y0(2)-1
            yinf(2)-0
            yinf(3)-0
            y0(5)-1
            yinf(5)-0
            y0(7)+m*y0(3)-0
            yinf(7)-0
            y0(10)+B_1*(1-y0(9))-0
            yinf(9)-0
            y0(12)+B_2*(1-y0(11))-0
            yinf(11)-0          
                ];
 end
 

The converted entropy generation number is

 NG = (1+(4/3)*Rd)*Re*y(10)*y(10)+2*(1+1/beta)*(2+K)*Br*Re*(1/omega_e)*(y(3)*y(3)+y(6)*y(6))+M*Br*Re*(1/omega_e)*(y(2)*y(2)+y(5)*y(5))+Re*Gamma*(omega_c/omega_e)*(omega_c/omega_e)*y(12)*y(12)+Re*Gamma*(omega_c/omega_e)*y(10)*y(12);

How to plot this equation within the code that I have.

[Nur Syahirah Wahid]

Hi. May I know if this has been solved?

[Saranya Shekar]

Hi, No this question has not been solved yet.