How to get 3 variables value from functions?

clc; clear all; close all; format short g; format compact;
global P A1 B1 C1 A2 B2 C2 Tc1 Pc1 Tc2 Pc2 R g12 g21 g11 g22 a12 x1 x2
P=1.01325*10^5; %Pa
A1=10.0331; B1=2940.46; C1=-35.93;
A2=9.5314; B2=2790.5; C2=-57.15;
Tc1=508.1; Pc1=47.01;
Tc2=523.2; Pc2=38.3;
R=8.314; %Pa.m^3/mol.K
g12=5; g21=5; g11=3; g22=2; a12=1.5;

x1=0:0.1:1;
x2=1-x1;

T=fzero(@NRTLFUN,300); %K

Error using fzero>localFirstFcnEval (line 729)
FZERO cannot continue because user-supplied function_handle ==> NRTLFUN failed with the error below.

Arrays have incompatible sizes for this operation.
Error in fzero (line 286)
fx = localFirstFcnEval(FunFcn,FunFcnIn,x,varargin{:});

plot(T,x1,y1)

function fT=NRTLFUN(T)
global P A1 B1 C1 A2 B2 C2 Tc1 Pc1 Tc2 Pc2 R g12 g21 g11 g22 a12 x1 x2
%Calculating lamda1 and lamda2
tau12=(g12-g22)/(R*T);
tau21=(g21-g11)/(R*T);
G12=exp(-a12*tau21);
G21=exp(-a12*tau12);
lamda1=exp((x2.^2).*((tau21.*((G21./(x1+(x2.*G21))).^2))+((G12.*tau12)./((x2+(x1.*G12)).^2))));
lamda2=exp((x1.^2).*((tau12.*((G12./(x2+(x1.*G12))).^2))+((G21.*tau21)./((x1+(x2.*G21)).^2))));

P1sat=exp(A1-(B1/(T+C1)))*10^5
P2sat=exp(A2-(B2/(T+C2)))*10^5

%fv = (P*v^3)-((R*T+P*b)*v^2)+(a*v)-(a*b);
a1 = (27*R^2*Tc1^2)/(64*Pc1)
b1 = (R*Tc1)/(8*Pc1)
a2 = (27*R^2*Tc2^2)/(64*Pc2)
b2 = (R*Tc2)/(8*Pc2)
fv1 = [P -(R*T+P*b1) a1 -a1*b1];
v1l = min(roots(fv1))
fv2 = [P -(R*T+P*b2) a2 -a2*b2];
v2l = min(roots(fv2))

y1=x1.*lamda1*P1sat*v1l/P
y2=x2.*lamda2*P2sat*v2l/P
fT=y1+y2-1;
end

How to solve this? I want to get T, y1 and y2 values. There is variable in range, which is x1. There are also some functions of T and one function include cubic function

[Rik]

Why are you using global variables? That makes debugging harder and is generally bad practice.

What have you tried yourself to find a solution? It looks like your function will error when you input a scalar. Have you confirmed it works as expected on scalar input?

Answer 1

[Sahil Jain]

Hi. From my understanding of the question, you want to write a function which outputs the values for "fT", "y1" and "y2". This can be done by changing the function declaration to:

function [fT,y1,y2]=NRTLFUN(T)

The function can then be called by:

[fT, y1, y2] = NRTLFUN(T);

However, this function will not work with "fzero" as "fzero" expects the input function to be scalar-valued. You can read more about function requirements for "fzero" on the documentation page.

[Nabilla Dewi Septiani]

Thank you! It is solved now