Search code examples
functionmatlabmatlab-figure

How can i use fsolve to plot the solutions to a function?


I have a variable of a that is equal to (weight./(1360*pi)).^(1/3), where the weight ranges between 4 and 8kg.

I then have guess of the time taken ,which is 14400 seconds.

The function in question is attached, where infinity is replaced by k=22. This function should be equal to 57/80

r/a can be replaced by 0.464, meaning that the multiplication of the summation can be written as 2/(0.464*pi).

alpha will be equal to 0.7*10^-7

How would i be able to plot the times taken for the masses to cook in hours, for weight in the given range?

I have tried to code this function for a couple of days now but it wont seem to work, due to array size issues and the general function just not working.


Solution

  • First, you need a master equation as a function of weight and t, which you want fsolve to find the zero of. Then for each weight, you can capture it in another function that you then solve for t:

    alpha = 0.7e-7;
    rbya = 0.464;
    k = 1:22;
    a = @(weight)(weight./(1360*pi)).^(1/3);
    eqn = @(weight,t)2/pi/rbya*sum((-1).^(k-1)./k.*sin(k*pi*rbya).*exp(-1.*k.^2.*pi^2.*alpha.*t./(a(weight).^2)))-57/80;
    
    weights = 4:8;
    ts = zeros(size(weights));
    for i = 1:numel(weights)
        sub_eqn = @(t)eqn(weights(i),t);
        ts(i)=fsolve(sub_eqn,14400);
    end
    
    plot(weights,ts/(60*60))
    xlabel("Weight (kg)")
    ylabel("Cooking Time (hrs)")
    

    If you want to solve the entire set of equations at once, then you need to be careful of array sizes (as you have experienced, read more here). k should be a column vector so that sum will sum along each column, and weights should be a row vector so that element-wise operations will repeat the k’s for each weight. You also need your list of initial guesses to be the same size as weights so that fsolve can have a guess for each weight:

    alpha = 0.7e-7;
    rbya = 0.464;
    k = (1:22)';
    
    a = @(weight)(weight./(1360*pi)).^(1/3);
    
    weights = 4:8;
    eqn = @(t)2/pi/rbya*sum((-1).^(k-1)./k.*sin(k*pi*rbya).*exp(-1.*k.^2.*pi^2.*alpha.*t./(a(weights).^2)))-57/80;
    
    ts=fsolve(eqn,repmat(14400,size(weights)));
    
    plot(weights,ts/(60*60))
    xlabel("Weight (kg)")
    ylabel("Cooking Time (hrs)")
    

    Note that you do get slightly different answers with the two methods.