diff --git a/C++/basicClasses/ngData_class.cc b/C++/basicClasses/ngData_class.cc index c75d6fbda3081308898a2435e13d1ca96624b273..8cc64ee70623c9a8617ef85cdbe0f61acb6487a2 100644 --- a/C++/basicClasses/ngData_class.cc +++ b/C++/basicClasses/ngData_class.cc @@ -606,6 +606,7 @@ void ngData_class<T, DIM_MESH>::Evaluate (const BaseMappedIntegrationRule & ir, double dval2; if (derivative == nullptr){ dval2 = this->Evaluate(ir[i]); + cout << "warning: no derivative set for dataClass" << endl; } else{ dval2 = derivative->Evaluate(ir[i]); @@ -613,10 +614,14 @@ void ngData_class<T, DIM_MESH>::Evaluate (const BaseMappedIntegrationRule & ir, // cout << "dval = " << dval << " =?= " << dval2 << endl; AutoDiff<1,double> res(val); - res.DValue(0) = dval2 ;//* x.DValue(0); + res.Value() = val; + res.DValue(0) = dval2;// * values(i, j).DValue(0); // cout << res.Value() << " " << res.DValue(0) << endl; values(i, j) = res; + + // cout << "i, j" << i << ", " << j << " " << ir[i].GetPoint()[0] << " " << values(i, j).Value() << " " << values(i, j).DValue(0) << endl; + } } } diff --git a/C++/nonLin/nonLinCF.cc b/C++/nonLin/nonLinCF.cc index cf936e825e64f42ee05112dc8248de550f93fd0f..2af32402a1a449dfba8ac000ab88a2c1caf09214 100644 --- a/C++/nonLin/nonLinCF.cc +++ b/C++/nonLin/nonLinCF.cc @@ -172,6 +172,12 @@ void nonLinCF<T, DIM>::Update(const shared_ptr<CoefficientFunction>& ptrInput){ } } else if(!std::strcmp(sFieldOut.c_str(), "mu_diff")){ + if (usePreisBiro && std::abs(val_in) <= concaveH){ + // double b = h*pastMu[el][p]; + // double hnew = bh->GetHfromB(b); + // return b / hnew; + val_in = oKL->EvaluateB(val_in * this->Get(el, ip).getxx()); + } this->Set(el, ip, oKL->DifferentiateH(val_in)); } else{ @@ -196,8 +202,16 @@ void nonLinCF<T, DIM>::Update(const shared_ptr<CoefficientFunction>& ptrInput){ this->Set(el, ip, oKL->EvaluateB(val_in)/val_in); } - else if(!std::strcmp(sFieldOut.c_str(), "nu_diff")) + else if(!std::strcmp(sFieldOut.c_str(), "nu_diff")){ + if (usePreisBiro && std::abs(val_in) > concaveB){ + // double h = b*pastMu[el][p][i]; + // b = bh->GetBfromH(h); + // return h / b; + + val_in = oKL->EvaluateH(val_in * this->Get(el, ip).getxx()); + } this->Set(el, ip, oKL->DifferentiateB(val_in)); + } else{ throw std::invalid_argument("invalid input parameter: valid parameters are (H, B, nu, nu_diff)"); } @@ -356,14 +370,14 @@ AutoDiff<1,double> nonLinCF<T, DIM>:: operator() (AutoDiff<1,double> x) const // double dval = (valr-vall) / (2*eps); if(strcmp(sFieldOut.c_str(), "B")){ - throw std::runtime_error("not implemented"); + throw std::runtime_error("to use NewtonMetho use H_KL, mu_KL and sFieldOut=B "); } double dval2 = this->oKL->DifferentiateH(x.Value()); // cout << "dval = " << dval << " =?= " << dval2 << endl; AutoDiff<1,double> res(val); - res.DValue(0) = dval2 ;//* x.DValue(0); + res.DValue(0) = dval2 * x.DValue(0); // cout << res.Value() << " " << res.DValue(0) << endl; diff --git a/C++/preisach/preisachScalar/everett.cc b/C++/preisach/preisachScalar/everett.cc index b24e920ced20c38c5e1f5a29e2101e20e6cf5c52..0f9370350c423d1d0e3dd9281e1ff49d0c60ea70 100644 --- a/C++/preisach/preisachScalar/everett.cc +++ b/C++/preisach/preisachScalar/everett.cc @@ -2311,6 +2311,8 @@ namespace everett_ns{ o << "isNonlin\t" << e.getIsNonLin() << endl; o << "isInverse\t" << e.getIsInverse() << endl; o << "dim\t" << e.getDIM() << endl; + // o << "dEFalling\t" << e.getDE_falling() << endl; + // o << "dERising\t" << e.getDE_rising() << endl; o << "evaluation version:\t"; switch (e.getEvalVersion()){ diff --git a/C++/preisach/preisachScalar/preisach.cc b/C++/preisach/preisachScalar/preisach.cc index 888010f596fc7dafdc338ca84ac91b36f494cb50..79b7fabb5641e7c8dc1c0c9d3add8688d2d7bdbb 100644 --- a/C++/preisach/preisachScalar/preisach.cc +++ b/C++/preisach/preisachScalar/preisach.cc @@ -1139,6 +1139,7 @@ namespace preisachScalar_ns{ #if VERSION_CALCDIFF==0 if(! conventionalMatDiff){ + // cout << "her" << endl; diff --git a/python/myPackage.py b/python/myPackage.py index 297031c514bcd338f724f7a40e9c147cd8e61b1c..4d3b0b9017e11178cb6158c875bc83a9f6f3029a 100644 --- a/python/myPackage.py +++ b/python/myPackage.py @@ -306,7 +306,7 @@ def myNonLinSolver(a, f, gfu, pre, eps=1e-3, Nmax=20, a_assemble=True, f_assembl if printrates: - print("it:", it, "stopcrit", stopcritval) + print("it:", it, "stopcrit", stopcritval, " > ", eps, " ? : ", stopcritval > eps) if stopcritval < eps or converge_func(gfu, gfu_o):