How to remove complex numbers from feval

8 views (last 30 days)
Ken Mercer
Ken Mercer on 19 Jan 2019
Commented: Ken Mercer on 20 Jan 2019
When I run this code to evaluate the equations (using simple feval functions) I get complex answers resulting in the plot shown. How does one handle this in the code to get a single line fit that I need and get rid of the complex solutions?
The 2 column table being read in is ejust xperimental data from uni thats been cleaned and Ive tried using statements with real(x) etc that doesnt work.
Many thanks for the help
GSI = 100;
D = 0;
Sig3 = TeoFQI.Sig3;
Sig1 = TeoFQI.Sig1;
HoekBrown = @(b,Sig3)(Sig3+b(1)*((b(2)*exp((GSI-100)/(28-14*D))*(Sig3/b(1))+ exp((GSI-100)/(9-3*D)))).^(0.5+1/6*(exp(-GSI/15)-exp(-20/3))));
opts = statset('nlinfit');
beta0 = [100;7];
beta = nlinfit(Sig3,Sig1,HoekBrown,beta0,opts)
S3min = min(Sig3)
S3max = max(Sig3)
figure()
plot(Sig3,Sig1,'o')
hold on
S3 = linspace(-6,S3max+4,100)
HB=@(Sig3x) Sig3x+beta(1)*((beta(2)*exp((GSI-100)./(28-14*D)).*(Sig3x/beta(1))+ exp((GSI-100)/(9-3*D)))).^(0.5+1/6.*(exp(-GSI/15)-exp(-20./3)));
y = feval(HB,S3)
plot(S3,real(y))
hold off
mb = beta(2).*exp((GSI-100)./(28-14*D));
s = exp((GSI-100)/(9-3*D));
a= 0.5+1/6*(exp(-GSI/15)-exp(-20/3));
SigciHold = beta(1)
SigN = @(Sig1x,Sig3x) ((Sig1x+Sig3x)/2 - (Sig1x-Sig3x)/2*(((1+a.*mb.*(mb*Sig3x./SigciHold+s).^(a-1))-1)/((1+a.*mb.*(mb*Sig3x./SigciHold+s).^(a-1))+1)));
SigN_Evaluated = feval(SigN,y,S3)
Tau = @(Sig1x,Sig3x) ((Sig1x-Sig3x).*(sqrt(1+a.*mb.*(mb.*Sig3./SigciHold+s).^(a-1)))./((1+a.*mb.*(mb.*Sig3./SigciHold+s).^(a-1))+1));
Tau_Evaluated = feval(Tau,y,S3)
figure()
plot(SigN_Evaluated,Tau_Evaluated,'o')
  4 Comments
Walter Roberson
Walter Roberson on 20 Jan 2019
nlinfit() does not appear to me to support bounds parameters.
It would help us to test if you supplied the data.
Ken Mercer
Ken Mercer on 20 Jan 2019
Walter many thanks you were corrent it is a bounds issue. Now resolved.

Sign in to comment.

Answers (0)

Categories

Find more on Graphics Performance in Help Center and File Exchange

Products


Release

R2018b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!