Joughin sliding-law inversion implemented and tested on the inverse e…

…xample. Seems to be working.
GHilmarG · Nov 28, 2023 · 0685e4a · 0685e4a
1 parent 413bf6f
commit 0685e4a
Show file tree

Hide file tree

Showing 8 changed files with 21 additions and 9 deletions.
diff --git a/BasalDrag.m b/BasalDrag.m
@@ -156,7 +156,7 @@
                 fprintf("Inversion using Coulomb sliding law not implemented. \n")
                 error("BasalDrag:InvalidCase","Inversion using Tsai sliding law not implemented.")
 
-            case {"Joughan","rCW-v0"}
+            case {"Joughin","rCW-v0"}
 
                  U=speed;
                  dFuvdC=   (U.^(1./m-1).*He.*(U+V0).^(-1./m).*(C+C0).^(-1./m-1))./m ;

diff --git a/CtrlVarValidityCheck.m b/CtrlVarValidityCheck.m
@@ -255,6 +255,15 @@
 
     end
 
+if CtrlVar.Inverse.MinimisationMethod=="MatlabOptimization-HessianBased" && CtrlVar.TriNodes>3
+
+
+    fprintf("Using CtrlVar.Inverse.MinimisationMethod=%s for other then linear elements (ie for CtrlVar.TriNodes>3) is not recommended. \n",CtrlVar.Inverse.MinimisationMethod)
+    fprintf('Consider setting CtrlVar.Inverse.MinimisationMethod="MatlabOptimization-GradientBased" or using linear elements. \n')
+    warning("UaInputs:ParameterCombinationNotRecommented","ParameterCombinationNotRecommented")
+
+end
+
 
 
     if contains(CtrlVar.Inverse.MinimisationMethod,'Hessian')
@@ -265,7 +274,7 @@
 
 end
 
-
+%%
 
 if isfield(CtrlVar,'AdaptMeshInterval')
 

diff --git a/CubicFit.m b/CubicFit.m
@@ -66,7 +66,7 @@
 
     if ab(1) == 0 || ~isreal(xmin) || isnan(xmin)
 
-        if ~isnan(Slope0)
+        if ~isnan(Slope)
             D=((y1-y0)/x1-Slope) ;
             xmin=-x1*Slope/2/D;
             if ~isnan(xmin)

diff --git a/GetInputsForInverseRun.m b/GetInputsForInverseRun.m
@@ -24,7 +24,7 @@
 [InvStartValues.AGlen,InvStartValues.n]=TestAGlenInputValues(CtrlVar,MUA,InvStartValues.AGlen,InvStartValues.n);
 %[Priors.AGlen,Priors.n]=TestAGlenInputValues(CtrlVar,MUA,Priors.AGlen,Priors.n);
 
-[InvStartValues.C,InvStartValues.m,InvStartValues.q,InvStartValues.muk]=TestSlipperinessInputValues(CtrlVar,MUA,InvStartValues.C,InvStartValues.m,InvStartValues.q,InvStartValues.muk);
+[InvStartValues.C,InvStartValues.m,InvStartValues.q,InvStartValues.muk]=TestSlipperinessInputValues(CtrlVar,MUA,InvStartValues.C,InvStartValues.m,InvStartValues.q,InvStartValues.muk,InvStartValues.V0);
 %[Priors.C,Priors.m]=TestSlipperinessInputValues(CtrlVar,MUA,Priors.C,Priors.m);
 
 %[Priors.rho,Priors.rhow]=TestDensityInputValues(CtrlVar,MUA,Priors.rho,Priors.rhow);

diff --git a/InversionValues.m b/InversionValues.m
@@ -8,6 +8,7 @@
         n=[];
         q=[];
         muk=[];
+        V0=[];
         C=[];
         AGlen=[];
         B=[];

diff --git a/PrintRunInfo.m b/PrintRunInfo.m
@@ -2,7 +2,7 @@
 function PrintRunInfo(CtrlVar)
 
 fprintf(CtrlVar.fidlog,'\n \n');
-fprintf(CtrlVar.fidlog,' **************************    %ca version: 8 November 2023 (alpha)  ********************** \n',218);
+fprintf(CtrlVar.fidlog,' **************************    %ca version: 28 November 2023 (alpha)  ********************** \n',218);
 fprintf(CtrlVar.fidlog,'    Run starts at %s  \n ',datetime);
 fprintf(CtrlVar.fidlog,'   Experiment: %s  \n ',CtrlVar.Experiment);
 

diff --git a/PriorProbabilityDistribution.m b/PriorProbabilityDistribution.m
@@ -16,6 +16,7 @@
         n=[];
         q=[];
         muk=[];
+        V0=[];
 
         AGlen=[];
         CovAGlen=[]

diff --git a/UaUtilities/PlotResultsFromInversion.m b/UaUtilities/PlotResultsFromInversion.m
@@ -859,16 +859,17 @@ function PlotResultsFromInversion(UserVar,CtrlVar,MUA,BCs,F,l,GF,InvStartValues,
             hold off
             yyaxis left
             semilogy(RunInfo.Inverse.Iterations,RunInfo.Inverse.J,'-bo','LineWidth',2)
-            ylabel('J','interpreter','latex')
+            ylabel('$J$','interpreter','latex')
 
             hold on
             semilogy(RunInfo.Inverse.Iterations,RunInfo.Inverse.I,'-gx')
-            ylabel('J and I')
+            ylabel('$J$ and $I$',Interpreter='latex')
 
             yyaxis right
             semilogy(RunInfo.Inverse.Iterations,RunInfo.Inverse.R,'-r+')
-            ylabel('R','interpreter','latex')
-            legend('Objective function','I','R','Location','southwest','interpreter','latex')
+            ylabel('$R$','interpreter','latex')
+            xlabel('Inverse iteration','interpreter','latex');
+            legend('Objective function','$I$','$R$','Location','southwest','interpreter','latex')
 
         end
-Original file line number
+Diff line change
@@ Expand Up / @@ -8,6 +8,7 @@ @@
             n=[];
             q=[];
             muk=[];
+            V0=[];
             C=[];
             AGlen=[];
             B=[];
@@ Expand Down @@