From 60f1f539970d8db89f08c3c504a81db9ddc80b7d Mon Sep 17 00:00:00 2001
From: Christophe Geuzaine <cgeuzaine@ulg.ac.be>
Date: Tue, 6 Mar 2018 16:33:48 +0100
Subject: [PATCH] sync with getdp
---
Team32/magstadyn_a.pro | 26 ++++++++++++++++++++------
Team32/param_EnergHyst.dat | 32 ++++++++++++++++++++++++++++----
Team32/param_EnergHyst.pro | 26 ++++++++++++++++++++++----
Team32/square_data.geo | 10 +++++++---
Team32/t32_data.geo | 9 +++++++--
5 files changed, 84 insertions(+), 19 deletions(-)
diff --git a/Team32/magstadyn_a.pro b/Team32/magstadyn_a.pro
index 6bad6b2..6aaeb80 100644
--- a/Team32/magstadyn_a.pro
+++ b/Team32/magstadyn_a.pro
@@ -70,6 +70,8 @@ Include "BH.pro"; // nonlinear BH caracteristic of magnetic material
Include "BH_anhysteretic.pro";
If(Flag_NL_law==NL_ENERGHYST)
Include "param_EnergHyst.pro";
+Else
+ N=0;
EndIf
@@ -209,12 +211,20 @@ Constraint {
{ Name MVP_2D ;
Case {
- { Region Corner ; Value 0.0 ; }
+ { Region Corner ; Value 0. ; }
{ Region SurfaceGe0 ; Type Assign ; Value 0. ; }
{ Region SurfaceGInf ; Type Assign ; Value 0. ; }
}
}
+ For i_ In {0:1}
+ { Name Dirichlet~{i_} ; Type Assign ;
+ Case {
+ { Region Boundary ; Type Assign; Value 0. * i_ ; }
+ }
+ }
+ EndFor
+
{ Name Current_2D ;
Case {
If(!Flag_CircuitCoupling)
@@ -255,7 +265,8 @@ FunctionSpace {
{ NameOfCoef ae1 ; EntityType NodesOf ; NameOfConstraint MVP_2D ; }
If (Flag_Degree_a == 2)
{ NameOfCoef ae2 ; // Only OK if homogeneous BC, otherwise specify zero-BC
- EntityType EdgesOf ; NameOfConstraint MagneticVectorPotential_2D ; }
+ //EntityType EdgesOf ; NameOfConstraint MagneticVectorPotential_2D ; }
+ EntityType EdgesOf ; NameOfConstraint Dirichlet_0 ; }
EndIf
}
}
@@ -433,7 +444,8 @@ Formulation {
// BF is constant per element (Vector + BF_Volume) => 1 integration point is enough
Galerkin { [ Dof{h} , {h} ] ;
In DomainNL ; Jacobian Vol ; Integration I1p ; }
- Galerkin { [ - (SetVariable[(h_Vinch_K[ {h}[1] ,
+ //Galerkin { [ -0*Vector[X[],Y[],Z[]] , {h} ] ; In DomainNL ; Jacobian Vol ; Integration I1p ; } //For Debug
+ Galerkin { [ - (SetVariable[(h_Vinch_K[ {h} , // new since 2/3/2018: use {h} instead of {h}[1] here (need to init bc by recomputing b from {h} in h_Vinch_K)
{d a}, {d a}[1],
{J_1}, {J_1}[1],
{J_2}, {J_2}[1],
@@ -487,7 +499,7 @@ Formulation {
// BF is constant per element (Vector + BF_Volume) => 1 integration point is enough
Galerkin { [ Dof{h} , {h} ] ;
In DomainNL ; Jacobian Vol ; Integration I1p ; }
- Galerkin { [ - (SetVariable[(h_Vinch_K[ {h}[1] ,
+ Galerkin { [ - (SetVariable[(h_Vinch_K[ {h} , // new since 2/3/2018: use {h} instead of {h}[1] here (need to init bc by recomputing b from {h} in h_Vinch_K)
{d a}, {d a}[1],
{J_1}, {J_1}[1],
{J_2}, {J_2}[1],
@@ -719,10 +731,12 @@ Resolution {
IterativeLoopN[ Nb_max_iter, RF_tuned[],
//*****Choose between one of the 3 following possibilities:*****
- //System { { A , Reltol_Mag, Abstol_Mag, Residual MeanL2Norm }} ]{ //1)
+ //System { { A , Reltol_Mag, Abstol_Mag, Solution MeanL2Norm }} ]{ //1a) //dx=x-xp; x=x
+ //System { { A , Reltol_Mag, Abstol_Mag, RecalcResidual MeanL2Norm }} ]{ //1b) //dx=res=b-Ax; x=b
+ //System { { A , Reltol_Mag, Abstol_Mag, Residual MeanL2Norm }} ]{ //1c) //dx=res=b-Ax; x=b #(default for square) #CHECK here
//PostOperation { { az_only , Reltol_Mag, Abstol_Mag, MeanL2Norm }} ]{ //2)
//PostOperation { { b_only , Reltol_Mag, Abstol_Mag, MeanL2Norm }} ]{ //3)
- PostOperation { { h_only , Reltol_Mag, Abstol_Mag, MeanL2Norm }} ]{ //4) //Need the above "INIT" for square with EnergHyst or JA because h_only = 0 at iter 1
+ PostOperation { { h_only , Reltol_Mag, Abstol_Mag, MeanL2Norm }} ]{ //4) //(default for t32) #CHECK here //Need the above "INIT" for square with EnergHyst or JA because h_only = 0 at iter 1
//**************************************************************
Evaluate[$res = $ResidualN, $resL = $Residual,$resN = $ResidualN,$iter = $Iteration-1];
diff --git a/Team32/param_EnergHyst.dat b/Team32/param_EnergHyst.dat
index e5e6931..cf500f0 100644
--- a/Team32/param_EnergHyst.dat
+++ b/Team32/param_EnergHyst.dat
@@ -1,7 +1,8 @@
+///*//default for t32 #CHECK here
dim00 = (Flag_3Dmodel==1)?3:2;
N00 = 3;
FLAG_TANORLANG00 = 1;
-FLAG_VARORDIFF00 = 1;
+FLAG_VARORDIFF00 = 2;
FLAG_INVMETHOD00 = 1;
FLAG_ROOTFINDING1D00 = 0;
FLAG_MINMETHOD00 = 1;
@@ -11,8 +12,31 @@ TOLERANCE_000 = 1.e-7;
TOLERANCE_NR00 = 1.e-7;
MAX_ITER_NR00 = 200;
TOLERANCE_OM00 = 1.e-11;
-MAX_ITER_OM00 = 700;
+MAX_ITER_OM00 = 1000;
FLAG_ANA00 = 1;
-TOLERANCE_NJ00 = 1.e-3;
+TOLERANCE_NJ00 = 1.e-14;
DELTA_000 = 1.e-0;
-FLAG_HOMO00 = 0;
\ No newline at end of file
+FLAG_HOMO00 = 1;
+//*/
+
+//default for square #CHECK here
+/*
+dim00 = (Flag_3Dmodel==1)?3:2;
+N00 = 3;
+FLAG_TANORLANG00 = 1;
+FLAG_VARORDIFF00 = 2;
+FLAG_INVMETHOD00 = 1;
+FLAG_ROOTFINDING1D00 = 0;
+FLAG_MINMETHOD00 = 1;
+FLAG_WARNING00 = 3 ;
+TOLERANCE_JS00 = 1.e-3;
+TOLERANCE_000 = 1.e-8;
+TOLERANCE_NR00 = 1.e-8;
+MAX_ITER_NR00 = 200;
+TOLERANCE_OM00 = 1.e-11;
+MAX_ITER_OM00 = 1000;
+FLAG_ANA00 = 1;
+TOLERANCE_NJ00 = 1.e-14;
+DELTA_000 = 1.e-0;
+FLAG_HOMO00 = 0;
+//*/
\ No newline at end of file
diff --git a/Team32/param_EnergHyst.pro b/Team32/param_EnergHyst.pro
index 80e4fe5..018a91a 100644
--- a/Team32/param_EnergHyst.pro
+++ b/Team32/param_EnergHyst.pro
@@ -82,7 +82,7 @@ FLAG_MINMETHOD = 6 --> steepest descent (gsl)
EndIf
If(FLAG_TANORLANG==2)
/*
- // Material M23535A
+ // Material M23535A (default for square) #CHECK here
Ja = 0.595;
ha = 4100;
Jb = 1.375;
@@ -100,7 +100,7 @@ FLAG_MINMETHOD = 6 --> steepest descent (gsl)
Jb=0.549518 ;
hb=188.516;
//*/
-///* //FH(xini=400) // NEW 1/9/2017
+///* //FH(xini=400) // NEW 1/9/2017 (default for t32) #CHECK here
Ja = 0.792234;
ha = 9.082095;
Jb = 0.790993;
@@ -115,7 +115,7 @@ hb = 399.474553851
EndIf
If(N==3)
///*
- // Material M25050A
+ // Material M25050A (default for square) #CHECK here
Js_1 = 0.11;
Js_2 = 0.8;
Js_3 = 0.31;
@@ -149,7 +149,7 @@ hb = 399.474553851
Js_2 = w_2*(Ja+Jb);
Js_3 = w_3*(Ja+Jb);
//*/
-/* //FH(xini=400) // NEW 1/9/2017 (chamonix)
+/* //FH(xini=400) // NEW 1/9/2017 (chamonix) (default for t32) #CHECK here
chi_1 = 0.;
chi_2 = 53.789872;
chi_3 = 233.910892;
@@ -338,6 +338,24 @@ Label_FLAG_ANA="...";
If (FLAG_MINMETHOD == 11 && FLAG_VARORDIFF==1)
Label_FLAG_MINMETHOD ="FLAG_MINMETHOD=11 --> NEW steepest descent naive (homemade)";
EndIf
+ If (FLAG_MINMETHOD == 22 && FLAG_VARORDIFF==1)
+ Label_FLAG_MINMETHOD ="FLAG_MINMETHOD=11 --> NEW conjugate Fletcher-Reeves (homemade)";
+ EndIf
+ If (FLAG_MINMETHOD == 33 && FLAG_VARORDIFF==1)
+ Label_FLAG_MINMETHOD ="FLAG_MINMETHOD=11 --> NEW conjugate Polak-Ribiere (homemade)";
+ EndIf
+ If (FLAG_MINMETHOD == 333 && FLAG_VARORDIFF==1)
+ Label_FLAG_MINMETHOD ="FLAG_MINMETHOD=11 --> NEW conjugate Polak-Ribiere+ (homemade)";
+ EndIf
+ If (FLAG_MINMETHOD == 1999 && FLAG_VARORDIFF==1)
+ Label_FLAG_MINMETHOD ="FLAG_MINMETHOD=11 --> NEW conjugate Dai Yuan 1999 (p.85) (homemade)";
+ EndIf
+ If (FLAG_MINMETHOD == 2005 && FLAG_VARORDIFF==1)
+ Label_FLAG_MINMETHOD ="FLAG_MINMETHOD=11 --> NEW conjugate Hager Zhang 2005 (p.161) (homemade)";
+ EndIf
+ If (FLAG_MINMETHOD == 77 && FLAG_VARORDIFF==1)
+ Label_FLAG_MINMETHOD ="FLAG_MINMETHOD=11 --> NEW newton (homemade)";
+ EndIf
If (FLAG_MINMETHOD == 2 && FLAG_VARORDIFF==1)
Label_FLAG_MINMETHOD ="FLAG_MINMETHOD=2 --> conjugate fr (gsl)";
EndIf
diff --git a/Team32/square_data.geo b/Team32/square_data.geo
index 8857a48..801c10f 100644
--- a/Team32/square_data.geo
+++ b/Team32/square_data.geo
@@ -49,7 +49,6 @@ Flag_NLRes00 = 0; // 0: use classical IterativeLoop to solve the NL (non linear)
// 1: use IterativeLoopN to solve the NL (non linear) problem
// 2: solve the NL (non linear) problem "by hand"
-
// non linear loop default parameters
Nb_max_iter00 = 50; // maximum number of NL (non linear) iterations
stop_criterion00 = 1e-5; // strop criterion for the NL (non linear) iteration
@@ -58,8 +57,13 @@ Flag_AdaptRelax00 = 1; // set to 1 to test different relaxation factors
relax_max00 = 1 ;
relax_min00 = 0.1;
relax_numtest00 = 10;
-TestAllFactors00 = 0; // 0 : stops when the residual goes up !!
- // 1 : try every relaxation factors and keep the optimal one
+TestAllFactors00 = 2; // 0 : try first relaxation factors (from the list) and stops when the residual goes up
+ // 1 : try every relaxation factors (from the list) and keep the optimal one
+ // 2 : find the maximum relaxation factor that decreases the residual:
+ // - start with the relaxation factor from the previous time step or from the previous iteration
+ // - the relaxation factor is multiplied by a ratio as long as the residual decreases
+ // - the relaxation factor is decreased by a ratio until a decreasing residual is found
+ // [3 : Build a parabola based on three relaxation factors and deduce a minimizing relaxation factor (NOT WORKING!!)]
Reltol_Mag00 = stop_criterion00; // 0 before with IterativeLoopN
Abstol_Mag00 = stop_criterion00;
Reltoldx_Mag00 = 1e-5*stop_criterion00;
diff --git a/Team32/t32_data.geo b/Team32/t32_data.geo
index 3a1c01c..bd2eb7d 100644
--- a/Team32/t32_data.geo
+++ b/Team32/t32_data.geo
@@ -44,8 +44,13 @@ Flag_AdaptRelax00 = 1; // set to 1 to test different relaxation factors
relax_max00 = 1 ;
relax_min00 = 0.1;
relax_numtest00 = 10;
-TestAllFactors00 = 1; // 0 : stops when the residual goes up !!
- // 1 : try every relaxation factors and keep the optimal one
+TestAllFactors00 = 1; // 0 : try first relaxation factors (from the list) and stops when the residual goes up
+ // 1 : try every relaxation factors (from the list) and keep the optimal one
+ // 2 : find the maximum relaxation factor that decreases the residual:
+ // - start with the relaxation factor from the previous time step or from the previous iteration
+ // - the relaxation factor is multiplied by a ratio as long as the residual decreases
+ // - the relaxation factor is decreased by a ratio until a decreasing residual is found
+ // [3 : Build a parabola based on three relaxation factors and deduce a minimizing relaxation factor (NOT WORKING!!)]
Reltol_Mag00 = stop_criterion00; // 0 before with IterativeLoopN
Abstol_Mag00 = stop_criterion00;
Reltoldx_Mag00 = 1e-5*stop_criterion00;
--
GitLab