gMuNuGeneric={{-f[x], 0 ,0,0}, {0, 1/f[x], 0, 0}, {0, 0, C[x]^2, 0}, {0, 0, 0,C[x]^2}};
MatrixForm[gMuNuGeneric];
GMuNuGeneric=Inverse[gMuNuGeneric];
coord2D={t, x, y, z};
sqrtngGeneric=Simplify[Sqrt[-Det[gMuNuGeneric]], {C[x]>0}];
box4DGeneric=0;box4D2Generic=0;
Do[ {box4D2Generic=1/sqrtngGeneric*(D[sqrtngGeneric*GMuNuGeneric[[i,j]]*D[Theta[t,x]/C[x],coord2D[[j]]] ,coord2D[[i]]]),box4DGeneric=box4DGeneric+box4D2Generic}, {i, 1, 2}, {j, 1, 2}];
box4DGenericDimple=Collect[Simplify[box4DGeneric], {Derivative[0,1][Theta][t,x], Derivative[2,0][Theta][t,x], Derivative[0,2][Theta][t,x],Theta[t,x]}];
box4DGenericDimpleCorrected=Simplify[box4DGenericDimple*C[x]];

The Potential contains already contains an overal negative sign because it is found through the Mode equation.

PotentialGeneric=Coefficient[box4DGenericDimpleCorrected,Theta[t,x]];
EffectivePotentialGeneric=PotentialGeneric*f[x];
A2[x_]=n[x]*(c[x]^2-v[x]^2)/c[x];
DA2[x_]=D[A2[x], x];
B2[x_]=(Sqrt[n[x]/c[x]]);
DB2[x_]=D[B2[x], x];
DDB2[x_]=D[DB2[x], x];
EffectivePotentialGeneric2019=Expand[(EffectivePotentialGeneric/.{f[x]->A2[x], f'[x]->DA2[x], C[x]->B2[x],C'[x]->DB2[x] ,C''[x]->DDB2[x] })];

rho1plus=264/100;
rho1minus=11;
XXX=27/10*10^(-7);
TTT=10^(-4);
Mplus=rho1minus/rho1plus;
uNISt=166053906660/100000000000*10^-27;(*1.66053906660*10^-27*);
hbarNIST=1054571817/1000000000*10^-34;(*1.054571817*10^-34*);
MRB87NISt=86909187/1000000*uNISt;(*86.909187*uNISt*);
CInt2016=25*10^-5;
CExt2016=57*10^-5;
XiInt=hbarNIST/(MRB87NISt*CInt2016);
XiExt=hbarNIST/(MRB87NISt*CExt2016);
Xi=Sqrt[XiInt*XiExt]

sigma=Sqrt[Mplus-1]/(Xi);
nZero[x1_]=(1/XXX)*rho1plus*(Mplus+(1-Mplus)/(Cosh[sigma*(xstart1)]^2));
vZero=1/(TTT*nZero[xstart1]);
cZero=(XXX/TTT)*Sqrt[(XXX*nZero[xstart1])/(rho1plus*rho1minus^2)];
xstartSolver=NSolve[vZero==cZero, xstart1, Reals, WorkingPrecision->125];

n2[x_]=(1/XXX)*rho1plus*(Mplus+(1-Mplus)/(Cosh[sigma*(x+xstart1/.xstartSolver[[1]])]^2));(*This contains approximately the location of c[x]=v[x]*)

dn2=D[n2[x], x]/n2[x];
ddn2=D[dn2, x]/n2[x];
v2=1/(TTT*n2[x]);
dv2=D[v2, x]/n2[x];
ddv2=D[dv2, x]/n2[x];
c2=(XXX/TTT)*Sqrt[(XXX*n2[x])/(rho1plus*rho1minus^2)];
dc2=D[c2, x]/n2[x];
ddc2=D[dc2,x]/n2[x];

N[(c2-v2)/.x->0]

N[D[(c2-v2), x]/.x->0]

NSolve[v2==c2, x, Reals]

xstartEnd1=N[-xstart1/.xstartSolver[[1]], 121]

(*ncvPlot=Plot[{ n2[x]*XXX*TTT, c2, v2}, {x, 0, xstartEnd1}, PlotLegends->{ "n[x]*X*T","c[x]", "v[x]"}, ImageSize->Large]*)

EffectivePotential2019Genericv1Values=HeavisideTheta[xstartEnd1-x]*EffectivePotentialGeneric2019/.{c[x]->c2, v[x]->v2,n[x]->n2[x],v'[x]->dv2,c'[x]->dc2,n'[x]->dn2,c''[x]->ddc2, v''[x]->ddv2  ,  n''[x]->ddn2};
(*EffectivePotential2019PlotInterior=Plot[-EffectivePotential2019Genericv1Values, {x, 0, 3/2*xstartEnd1}, ImageSize->Large, ImageSize->Large, PlotRange->All, PlotLegends->{ "V_Eff"}]*)

(*Export["/Users/richarddudley/Dropbox/BH/BEC/nofxcofxvofx/Writeup/3/EffectivePotential2019PlotInterior.jpg",EffectivePotential2019PlotInterior]*)

EffectivePotential2019Genericv1Values2=EffectivePotentialGeneric2019/.{c[x]->c2, v[x]->v2,n[x]->n2[x],v'[x]->dv2,c'[x]->dc2,n'[x]->dn2,c''[x]->ddc2, v''[x]->ddv2  ,  n''[x]->ddn2};

N[c2/.x->-5000*Xi]

N[v2/.x->-5000*Xi]

N[c2/.x->xstartEnd1]

N[(1/(2*v2)*D[c2^2-v2^2, x])/.x->0, 50]

N[(D[c2-v2, x])/.x->0]

Extx1=-10^-22;
Extx3=Extx1;
Extx2=-Extx1;
Extx4=Extx2;
Exth=FullSimplify[-c2/(c2^2-v2^2)];
ExtDh=FullSimplify[D[Exth, x]];
Extxstar[y_]:=NIntegrate[(Exth), {x,  Extx3,y}, WorkingPrecision->80, PrecisionGoal->50, MaxRecursion->50]
Extxinitial1=N[-20*10^-6, 70];
Extxfinal=N[-10^-22,70];
Extxfinal2=N[-10^-22,70];
Extxstar[Extxinitial1]

x2=10^-22;
x4=x2;
h=FullSimplify[-c2/(c2^2-v2^2)];
Dh=FullSimplify[D[h, x]];
dxstardxAtxValue=N[(h)/. x->xstartEnd1, 80];
dxstarvalue1=NIntegrate[((h)), {x,  x2,xstartEnd1},MaxRecursion->1000,  WorkingPrecision->80, PrecisionGoal->50, MaxRecursion->100];
xstar[y_]:=(HeavisideTheta[xstartEnd1-y]*(NIntegrate[((h)), {x,  x2,y},MaxRecursion->1000,  WorkingPrecision->80, PrecisionGoal->50, MaxRecursion->100])+N[HeavisideTheta[y-xstartEnd1]*((y-xstartEnd1)*dxstardxAtxValue+dxstarvalue1), 120]);

xstarSimple[y_]:=N[((y-xstartEnd1)*dxstardxAtxValue+dxstarvalue1), 80];

ExteriorxstarValue=N[Extxstar[Extxinitial1], 80]

xstar[20*Xi]

xstarSimple[65*Xi]

xinitial1=10^-22;

xfinal=N[xstartEnd1,80];
xfinal2=N[xstartEnd1,80];
xstarFinal1=N[xstarSimple[xstartEnd1], 121];
xstarFinal2=N[xstarSimple[xstartEnd1], 121]
Vcheck1=N[-EffectivePotential2019Genericv1Values/.x->xinitial1, 60]
Vcheck2=N[-EffectivePotential2019Genericv1Values/.x->xfinal]
Vcheck2=N[-EffectivePotential2019Genericv1Values/.x->xfinal2]

xstartEnd1

xfinal

xstarFinal2-ExteriorxstarValue

xinitial has to be the same as x2  for this to be correct

TIRIofOmegaVeff=Table[{(10^(i/11)*10^(-6)),
s1Veff=NDSolve[{1/(h^2)*y1Veff''[x]-(Dh)/(h^3)*y1Veff'[x]-(-((10^(i/11)*10^(-6)))^2-EffectivePotential2019Genericv1Values2)*y1Veff[x]==0,y1Veff[xinitial1]==(Exp[I*(10^(i/11)*10^(-6))*xstar[xinitial1]]),y1Veff'[xinitial1]==N[((I*(10^(i/11)*10^(-6))*(h))/.x->xinitial1)*(Exp[I*(10^(i/11)*10^(-6))*xstar[xinitial1]]), 121]},y1Veff,{x,xinitial1,xstartEnd1}, WorkingPrecision->70, PrecisionGoal->20, AccuracyGoal->20, MaxSteps->10000000,InterpolationOrder->All, Method->"ExplicitRungeKutta"];

N[ (Evaluate[Exp[-I*(10^(i/11)*10^(-6))*(xstarFinal2-ExteriorxstarValue)]]/2*(Evaluate[y1Veff[xfinal2]/.s1Veff]-(I/((10^(i/11)*10^(-6))*dxstardxAtxValue)/.x->xfinal2)*Evaluate[y1Veff'[xfinal2]/.s1Veff])), 30],
N[(Evaluate[Exp[I*(10^(i/11)*10^(-6))*(xstarFinal2-ExteriorxstarValue)]]/2*(Evaluate[y1Veff[xfinal2]/.s1Veff]+(I/((10^(i/11)*10^(-6))*dxstardxAtxValue)/.x->xfinal2)*Evaluate[y1Veff'[xfinal2]/.s1Veff])), 30]}, {i,1,50}];

N[ (Evaluate[Exp[I*(10^(i/11)*10^(-6))*(xstarFinal2-ExteriorxstarValue)]]/2*(Evaluate[y1Veff[xfinal2]/.s1Veff]+(I/((10^(i/11)*10^(-6))*dxstardxAtxValue)/.x->xfinal2)*Evaluate[y1Veff'[xfinal2]/.s1Veff])), 30]

Evaluate[y1Veff'[xfinal2]/.s1Veff]

N[(y1Veff[xstartEnd1]/.s1Veff), 60]

xstar[xfinal2]

(((xfinal2-xstartEnd1)*dxstardxAtxValue+dxstarvalue1))

(*-EffectivePotential2019Genericv1Values*)

Omega1Veff=Flatten[TIRIofOmegaVeff[[All, 1]]];
TI1Veff=Flatten[Map[ Abs,TIRIofOmegaVeff[[All, 2]]]];
RI1Veff=Flatten[Map[Abs, TIRIofOmegaVeff[[All, 3]]]];
dataTIVeff=Transpose[{Omega1Veff, TI1Veff}];
dataRIVeff=Transpose[{Omega1Veff, RI1Veff}];

AlphaPlotVeff=ListLogLinearPlot[Re[dataTIVeff],PlotStyle->Automatic,AxesLabel->{"ω"}, PlotLegends->{"|α|"}, PlotRange->All, ImageSize->Large]

Export["E:\\Dropbox\\BH\\BEC\\nofxcofxvofx\\CorrectCVN\\2020Files\\Potentialvariations\\Plots\\M4p17AlphaPlotVeff50p1.jpg",AlphaPlotVeff,ImageResolution->200]

BetaPlotVeff=ListLogLinearPlot[Re[dataRIVeff],PlotStyle->Automatic,AxesLabel->{"ω"}, PlotLegends->{"|β|"}, PlotRange->All, ImageSize->Large]

Export["D:\\Dropbox\\BH\\BEC\\nofxcofxvofx\\CorrectCVN\\2020Files\\Plots\\M4p17BetaPlotVeff50p1.jpg",BetaPlotVeff,ImageResolution->200]

Abs[dataTIVeff[[ 1, 2]]]^2-Abs[dataRIVeff[[ 1, 2]]]^2

Abs[dataTIVeff[[ 2, 2]]]^2-Abs[dataRIVeff[[ 2, 2]]]^2

Abs[dataTIVeff[[ 24, 2]]]^2-Abs[dataRIVeff[[ 19, 2]]]^2

Plot[1/(HeavisideTheta[xstartEnd1-x]*h^2+HeavisideTheta[x-xstartEnd1]*dxstardxAtxValue^2), {x, 10^-2, 1}]

Plot[(HeavisideTheta[xstartEnd1-x]*Dh)/(HeavisideTheta[xstartEnd1-x]*h^3+HeavisideTheta[x-xstartEnd1]*dxstardxAtxValue^3), {x, 10^-2, 1}, PlotRange->All]

Export["E:\\Dropbox\\BH\\BEC\\nofxcofxvofx\\CorrectCVN\\2020Files\\Potentialvariations\\RawData\\M4p17alphabetaTableData50p1.mx",TIRIofOmegaVeff]