Plotting a Piecewise function Dg?

Hello community I look to get the figure of Dg, please verifiy if the definition of Piecewise function Dg is correct and Thank you
clc
clear all
x=[]; y=[]; z=[];
for n=1:1001
x1=0.01*(n-1);
x2=4;
B=[0,x1,4,x2;
x1,3,x2,x2;
x1,0,5,x2;
0,5,x2,x1];
%{
Pd=eig(B);
if max(real(Pd))<0
disp('fail');
disp(n);
end
%}
B11=[B(1,1),B(1,2);B(2,1),B(2,2)];
B22=[B(3,3),B(3,4);B(4,3),B(4,4)];
B12=[B(1,3),B(1,4);B(2,3),B(2,4)];
V1_2=[B11,B12;B12',B22];
A1=det(B11);
B1=det(B22);
C1=det(B12);
D1=det(V1_2);
Sum=det(B11)+det(B22)+2.*det(B12);
nup=sqrt(Sum+sqrt(Sum.^2-4.*det(V1_2)))./sqrt(2);
nun=sqrt(Sum-sqrt(Sum.^2-4.*det(V1_2)))./sqrt(2);
%defin quantum discord
fmp=(nup+1)/2*log((nup+1)/2)-(nup-1)/2*log((nup-1)/2);
fmn=(nun+1)/2*log((nun+1)/2)-(nun-1)/2*log((nun-1)/2);
fB1=(sqrt(B1)+1)/2*log((sqrt(B1)+1)/2)-(sqrt(B1)-1)/2*log((sqrt(B1)-1)/2);
Dg1=fB1-fmp-fmn;
%inf epsilon
G1=(2*C1^2+(B1-1)*(D1-A1)+2*abs(C1)*sqrt(C1^2+(-1+B1)*(-A1+D1)))/(B1-1)^2;%val1
G2=(A1*B1-C1^2+D1-sqrt(C1^4+(-A1*B1+D1)^2-2*C1*(A1*B1+D1)))/2*B1;%val2
G3=(D1-A1*B1)^2-(1+B1)*C1^2*(A1+D1);%<=0 condition
DgG1=(sqrt(G1)+1)/2*log((sqrt(G1)+1)/2)-(sqrt(G1)-1)/2*log((sqrt(G1)-1)/2); %val if %G3<=0 condition
DgG2=(sqrt(G2)+1)/2*log((sqrt(G2)+1)/2)-(sqrt(G2)-1)/2*log((sqrt(G2)-1)/2); %else G3>0
if G3 <= 0
Dg = real(DgG1+Dg1);
else
Dg = real(DgG2+Dg1);
end
x(n)=x1; y(n)=Dg;
n=n+1;
end
plot(x,y)

3 Comments

How can we know if the function is correct, sicne you never tell us what the function should be in terms of mathematics?
Of course it is perfectly correct. It does EXACTLY what the code tells it to do. If it should work differently, then you should have written different code.
I need to define this function which is given in the code as G1 and G2 and the condition as G3. I want to plot this function as a function of other parameters such as temperature and coupling, since the variables A, B, C, and D are indirectly dependent on the mentioned parameters.@John D'Errico

Sign in to comment.

Answers (0)

Categories

Find more on Mathematics in Help Center and File Exchange

Products

Release

R2019a

Tags

Asked:

on 28 Apr 2023

Commented:

on 28 Apr 2023

Community Treasure Hunt

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

Start Hunting!