-2

enter image description here

% Condenser and hated fluid Data
UAc = 30600 ;   %W/K
mc = 6.8 ;      %kg/s

% Evaporator Data
UAe = 26500 ;   %W/K
me = 7.6 ;      %kg/s
%input data
ta = input('Input Temperature at evaporator inlet water(celsius),ta : ');
tb = input('Input Temperature at condensor  inlet water(celsius),tb : ');
cp = input('Specific heat capacity (kJ/kgK),cp : ');
tol = 0.0001 %tol is converge tolerance

%guess
te= 5;
tc= 45;
qe= 160;
qc= 190;
P = 40;
iter = 0 ; %iteration
xnew = [te;tc;qe;qc;P];
xold = zeros(size(xnew));

while norm( xnew- xold) > tol
    iter = iter+1;
    xold = xnew ;
    %update te,tc,qe,qc,P
    te= xold(1);
    tc= xold(2);
    qe= xold(3);
    qc= xold(4);
    P = xold(5);

    %Create equations for differrential;

    %Refigeration capacity
    f1 = sym('239.5+(10.073*te)-(0.109*(te^2))-(3.41*tc)-(0.0025*(tc^2))-(0.203*te*tc)+(0.0082*(te^2)*tc)+(0.0013*te*(tc^2))-(0.000080005*(te^2)*(tc^2))- qe');

    %Power required by compressor
    f2 = sym('-2.634-(0.3081*te)-(0.00301*(te^2))+(1.066*tc)-(0.00528*(tc^2))-(0.0011*te*tc)-(0.000306*(te^2)*tc)+(0.000567*te*(tc^2))-(0.0000031*(te^2)*(tc^2))-P');

    %Evaporator heat transfer rate
    f3 = sym('me*cp*(ta-te)*(1-exp(UAe/(me*cp*1000)))- qe');

    %Condensor heat transfer rate
    f4 = sym('mc*cp*(tc-tb)*(1-exp(UAe/(mc*cp*1000)))- qc');

    %Energy balance
    f5 = sym('P+qe-qc');

    % Newton raphson method

    %diff f1
    f1te = diff(f1,'te'); % differential f1 by te
    f1tc = diff(f1,'tc'); % differential f1 by tc
    f1qe = diff(f1,'qe'); % differential f1 by qe
    f1qc = diff(f1,'qc'); % differential f1 by qc
    f1P = diff(f1,'P'); % differential f1 by P

    %diff f2
    f2te = diff(f2,'te'); % differential f2 by te
    f2tc = diff(f2,'tc'); % differential f2 by tc
    f2qe = diff(f2,'qe'); % differential f2 by qe
    f2qc = diff(f2,'qc'); % differential f2 by qc
    f2P = diff(f2,'P'); % differential f2 by P

    %diff f3
    f3te = diff(f3,'te'); % differential f3 by te
    f3tc = diff(f3,'tc'); % differential f3 by tc
    f3qe = diff(f3,'qe'); % differential f3 by qe
    f3qc = diff(f3,'qc'); % differential f3 by qc
    f3P = diff(f3,'P'); % differential f3 by P

    %diff f4
    f4te = diff(f4,'te'); % differential f4 by te
    f4tc = diff(f4,'tc'); % differential f4 by tc
    f4qe = diff(f4,'qe'); % differential f4 by qe
    f4qc = diff(f4,'qc'); % differential f4 by qc
    f4P = diff(f4,'P'); % differential f4 by P

    %diff f5
    f5te = diff(f5,'te'); % differential f5 by te
    f5tc = diff(f5,'tc'); % differential f5 by tc
    f5qe = diff(f5,'qe'); % differential f5 by qe
    f5qc = diff(f5,'qc'); % differential f5 by qc
    f5P = diff(f5,'P'); % differential f5 by P

    J = [f1te f1tc f1qe f1qc f1P;f2te f2tc f2qe f2qc f2P;f3te f3tc f3qe f3qc f3P;f4te f4tc f4qe f4qc f4P;f5te f5tc f5qe f5qc f5P];

    xnew = xold - J\[te;tc;qe;qc;P];
    disp(fprintf('iter=%6.5f,  te=%6.5f,  tc=%6.5f, qe=%6.5f, qc=%6.5f, P=%6.5f', iter,xnew));

end
rozsasarpi
  • 1,592
  • 18
  • 33

1 Answers1

1

The error message seems to be clear. Why don't you just display the symbolic object xnew without fprint:

disp(fprintf('iter=%6.5f,  te=%6.5f,  tc=%6.5f, qe=%6.5f, qc=%6.5f, P=%6.5f', iter));
disp(xnew)

On the other hand based on your code I assume xnew should be a numeric vector instead of a symbolic variable. This is because your J variable is symbolic. I do not know what you want to accomplish because no description or even a question was given.. but you can try to substitute (subs) into J to get a numeric vector. It is not a good idea to use the same variable name in symbolic expressions and for numeric data, e.g. tc, tn...

You should not hate fluids.

rozsasarpi
  • 1,592
  • 18
  • 33