implementation of the bisection method

Hi everyone,

I'd like to use the bisection method to find the V values that will make my set of functions equal to zero (V_ls).

func=@(V) -Pstar+(R*T*(((1-c)/V)-((c/b)*log(1-b/V)))-a*c/V^2)

note that a,b,c,R retain always the same value, the difference that makes my set of equations is the values of Pstar,T that will then determine the V value that I'm looking for.

I successfully managed to do this process for one given function, by setting the value of Pstar and T.

I'd like to be able to repeat this process automatically for different values of Pstar and T that are stored in a vector.

Pstar=Pstars(i)
T=temp(i)
%i have 30 values (given by the problem) in each vector, i.e. i=30

%range values where the solution can be found.
V=0.0320;
lim_max= 0.1951;

%bisection method

number_iterations= 100 % i can set this value arbitrarily
for i=1:number_iterations
    Vls=(V+lim_max)/2;
    if fun(Vls)<0
        
        lim_max=Vls;
    else
        V=Vls;
    end
end

I suppose I should put this for cicle into a new one, but I really cant manage to get the results I want:

an output vector V_ls that contains each solution given by the bisection method for each of my 30 equations.

Many thanks !

[Torsten]

I don't see any executable code for one pair (T,Pstar).

[Ivan De Fazio]

thanks for your answer !

English isnt my first language so I think I will explain a little bit further just in case I misunderstood your answer.

I know this process works for all 30 equations

what I wish is for the code to iterate the same process for all 30 of them.

for example, if I am at iteration 2, I want to automatically input the value of temp(2) and Pstar(2) (i know both numerical values since they’re given by the problem) into my equation, then i’d like to run the bisection method and get the corresponding Vls and so forth. at the end of the whole cycle, when all 30 equations have been processed I want to generate a vector that contains all 30 Vls

[Torsten]

Yes, I understand this. But the code above does not work for one pair. How should I change it for 30 pairs then ?

If you want to try it on your own:

Make a function "bisection" that is called with the necessary inputs (temp(i),Pstar(i),fun,maximum number of iterations,...) and returns V(i). Then call this function in a loop like:

for i=1:30
    Vsol(i) = bisection(temp(i),Pstar(i),fun,Vinit,max_iter,...);
end

[Ivan De Fazio]

thanks for your help and patience.

This is my code that works for one iteration:

clear all;
close all;
clc;

%fixed inputs
c=25.8972;
b=0.0320;
a=0.7535;
R=0.0821;
Vc=0.1951;

%chosen couple T-Pstar
T=385.15;
Pstar=36.6133;

fun=@(V) -Pstar+(R*T*(((1-c)/V)-((c/b)*log(1-b/V)))-a*c/V^2);
n=100;

V=b;
limmax=Vc;

for i=1:n
    Vls=(V+limmax)/2;
    if fun(Vls)<0
        limmax=Vls;
    else 
        V=Vls;    
    end 
end

sol=Vls;
disp(sol)

    0.1038

i've tried for over an hour to create a function just as you said but it keeps giving me the same result.

I've been stuck with this for days.

thanks again

Answer 1

[Torsten]

%fixed inputs
c=25.8972;
b=0.0320;
a=0.7535;
R=0.0821;
Vc=0.1951;
% Generate 30 artificial (Pstar,T) pairs to test
Pstar_array = 36.6133*(1.05).^(1:30);
T_array = 385.15+(1:30);
% Call bisection for 30 (Pstar,T) pairs
for i=1:numel(T_array)
    Pstar = Pstar_array(i);
    T = T_array(i);
    fun = @(V) -Pstar+(R*T*(((1-c)/V)-((c/b)*log(1-b/V)))-a*c/V^2);
    V(i) = bisection(fun,b,Vc);
    error(i) = fun(V(i));
end
figure(1)
plot(1:30,V)

figure(2)
plot(1:30,error)

function sol = bisection(fun,V0,Vc)
n=100;
V=V0;
limmax=Vc;
for i=1:n
    Vls=(V+limmax)/2;
    if fun(Vls)<0
        limmax=Vls;
    else 
        V=Vls;    
    end 
end
sol=Vls;
end