1 #ifndef ROL_DIODECIRCUIT_HPP
2 #define ROL_DIODECIRCUIT_HPP
38 Teuchos::RCP<std::vector<Real> >
Imeas_;
40 Teuchos::RCP<std::vector<Real> >
Vsrc_;
59 Objective_DiodeCircuit(Real Vth, Real Vsrc_min, Real Vsrc_max, Real Vsrc_step, Real true_Is, Real true_Rs,
bool lambertw, Real noise,
bool use_adjoint,
int use_hessvec){
64 int n = (Vsrc_max-Vsrc_min)/Vsrc_step + 1;
65 Vsrc_ = Teuchos::rcp(
new std::vector<Real>(n,0.0));
66 Imeas_ = Teuchos::rcp(
new std::vector<Real>(n,0.0));
67 std::ofstream output (
"Measurements.dat");
68 Real left = 0.0; Real right = 1.0;
71 std::cout <<
"Generating data using Lambert-W function." << std::endl;
73 (*Vsrc_)[i] = Vsrc_min+i*Vsrc_step;
76 (*Imeas_)[i] += noise*pow(10,
int(log10((*
Imeas_)[i])))*(( (Real)rand() / (Real)RAND_MAX ) * (right - left) + left);
80 output << std::setprecision(8) << std::scientific << (*Vsrc_)[i] <<
" " << (*Imeas_)[i] <<
"\n";
86 std::cout <<
"Generating data using Newton's method." << std::endl;
88 (*Vsrc_)[i] = Vsrc_min+i*Vsrc_step;
90 (*Imeas_)[i] =
Newton(I0,Vsrc_min+i*Vsrc_step,true_Is,true_Rs);
92 (*Imeas_)[i] += noise*pow(10,
int(log10((*
Imeas_)[i])))*(( (Real)rand() / (Real)RAND_MAX ) * (right - left) + left);
96 output << std::setprecision(8) << std::scientific << (*Vsrc_)[i] <<
" " << (*Imeas_)[i] <<
"\n";
118 for(
int k=0;std::getline(input_file,line);++k){dim=k;}
120 input_file.seekg(0,std::ios::beg);
121 Vsrc_ = Teuchos::rcp(
new std::vector<Real>(dim,0.0));
122 Imeas_ = Teuchos::rcp(
new std::vector<Real>(dim,0.0));
124 std::cout <<
"Using input file to generate data." <<
"\n";
125 for(
int i=0;i<dim;i++){
129 (*Imeas_)[i] = Imeas;
149 Real
diode(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
150 return I-Is*(exp((Vsrc-I*Rs)/
Vth_)-1);
154 Real
diodeI(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
155 return 1+Is*exp((Vsrc-I*Rs)/
Vth_)*(Rs/
Vth_);
159 Real
diodeIs(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
160 return 1-exp((Vsrc-I*Rs)/
Vth_);
164 Real
diodeRs(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
165 return Is*exp((Vsrc-I*Rs)/
Vth_)*(I/
Vth_);
169 Real
diodeII(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
174 Real
diodeIIs(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
175 return exp((Vsrc-I*Rs)/
Vth_)*(Rs/
Vth_);
179 Real
diodeIRs(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
184 Real
diodeIsIs(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
189 Real
diodeIsRs(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
190 return exp((Vsrc-I*Rs)/
Vth_)*(I/
Vth_);
194 Real
diodeRsRs(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
205 Real
Newton(
const Real I,
const Real Vsrc,
const Real Is,
const Real Rs){
210 Real fval =
diode(IN,Vsrc,Is,Rs);
215 for (
int i = 0; i < MAXIT; i++ ) {
216 if ( std::fabs(fval) < TOL ) {
220 dfval =
diodeI(IN,Vsrc,Is,Rs);
221 if( std::fabs(dfval) < EPS ){
222 std::cout <<
"denominator is too small" << std::endl;
227 IN_tmp = IN - alpha*fval/dfval;
228 fval_tmp =
diode(IN_tmp,Vsrc,Is,Rs);
229 while ( std::fabs(fval_tmp) >= (1.0-1.e-4*alpha)*std::fabs(fval) ) {
231 IN_tmp = IN - alpha*fval/dfval;
232 fval_tmp =
diode(IN_tmp,Vsrc,Is,Rs);
233 if ( alpha < std::sqrt(EPS) ) {
278 void lambertw(Real x, Real &w,
int &ierr, Real &xi){
280 const Real turnpt = -exp(-1.), c1 = 1.5, c2 = .75;
281 Real r, r2, r3, s, mach_eps, relerr = 1., diff;
292 if( x == 0. )
return;
293 if( x < (1-c2) ) w = x*(1.-x + c1*x*x);
299 w = x*(1.0-x + c1*x*x);
300 xi = log(1.0-x + c1*x*x) - w;
304 if( diff < 0.0 ) diff = -diff;
305 w = -1 + sqrt(2.0*exp(1.))*sqrt(x-turnpt);
306 if( diff == 0.0 )
return;
318 while( relerr > mach_eps && i<maxit){
322 s = 6.*(w+1.0)*(w+1.0);
323 w = w * ( 1.0 + r + r2/(2.0*( w+1.0)) - (2. * w -1.0)*r3/s );
324 if( w * x < 0.0 ) w = -w;
333 if(relerr < 0.0 ) relerr = -relerr;
336 if( i == maxit ) ierr = 2;
352 Real arg1 = (Vsrc + Is*Rs)/
Vth_;
353 Real evd = exp(arg1);
354 Real lambWArg = Is*Rs*evd/
Vth_;
358 lambertw(lambWArg, lambWReturn, ierr, lambWError);
359 if(ierr == 1){std::cout <<
"LambertW error: argument is not in the domain" << std::endl;
return -1.0;}
360 if(ierr == 2){std::cout <<
"LambertW error: BUG!" << std::endl;}
361 Real Id = -Is+
Vth_*(lambWReturn)/Rs;
369 Teuchos::RCP<std::vector<Real> > Ip =
370 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(I)).getVector());
371 Teuchos::RCP<const std::vector<Real> > Sp =
379 for(
int i=0;i<n;i++){
387 for(
int i=0;i<n;i++){
388 (*Ip)[i] =
Newton(I0,(*
Vsrc_)[i],(*Sp)[0],(*Sp)[1]);
402 Teuchos::RCP<const std::vector<Real> > Sp =
406 Teuchos::RCP<std::vector<Real> > Ip =
407 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(I)).getVector());
413 for(
int i=0;i<n;i++){
414 val += ((*Ip)[i]-(*Imeas_)[i])*((*Ip)[i]-(*Imeas_)[i]);
428 Teuchos::RCP<std::vector<Real> > lambdap =
429 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(lambda)).getVector());
431 Teuchos::RCP<const std::vector<Real> > Ip =
433 Teuchos::RCP<const std::vector<Real> > Sp =
437 for(
int i=0;i<n;i++){
438 (*lambdap)[i] = ((*Imeas_)[i]-(*Ip)[i])/
diodeI((*Ip)[i],(*Vsrc_)[i],(*Sp)[0],(*Sp)[1]);
451 Teuchos::RCP<std::vector<Real> > sensp =
452 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(sens)).getVector());
453 Teuchos::RCP<const std::vector<Real> > Ip =
455 Teuchos::RCP<const std::vector<Real> > Sp =
459 for(
int i=0;i<n;i++){
473 Teuchos::RCP<std::vector<Real> > sensp =
474 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(sens)).getVector());
475 Teuchos::RCP<const std::vector<Real> > Ip =
477 Teuchos::RCP<const std::vector<Real> > Sp =
481 for(
int i=0;i<n;i++){
489 Teuchos::RCP<std::vector<Real> > gp =
490 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(g)).getVector());
491 Teuchos::RCP<const std::vector<Real> > Sp =
497 Teuchos::RCP<std::vector<Real> > Ip =
498 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(I)).getVector());
505 StdVector<Real> lambda( Teuchos::rcp(
new std::vector<Real>(n,0.0) ) );
506 Teuchos::RCP<std::vector<Real> > lambdap =
507 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(lambda)).getVector());
513 (*gp)[0] = 0.0; (*gp)[1] = 0.0;
514 for(
int i=0;i<n;i++){
515 (*gp)[0] +=
diodeIs((*Ip)[i],(*
Vsrc_)[i],(*Sp)[0],(*Sp)[1])*(*lambdap)[i];
516 (*gp)[1] +=
diodeRs((*Ip)[i],(*
Vsrc_)[i],(*Sp)[0],(*Sp)[1])*(*lambdap)[i];
521 StdVector<Real> sensIs( Teuchos::rcp(
new std::vector<Real>(n,0.0) ) );
522 StdVector<Real> sensRs( Teuchos::rcp(
new std::vector<Real>(n,0.0) ) );
527 Teuchos::RCP<std::vector<Real> > sensIsp = Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(sensIs)).getVector());
528 Teuchos::RCP<std::vector<Real> > sensRsp = Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(sensRs)).getVector());
531 std::ofstream output (
"Sensitivities.dat");
532 for(
int k=0;k<n;k++){
533 if(output.is_open()){
534 output << std::scientific << (*sensIsp)[k] <<
" " << (*sensRsp)[k] <<
"\n";
539 (*gp)[0] = 0.0; (*gp)[1] = 0.0;
540 for(
int i=0;i<n;i++){
541 (*gp)[0] += ((*Ip)[i]-(*Imeas_)[i])*(*sensIsp)[i];
542 (*gp)[1] += ((*Ip)[i]-(*Imeas_)[i])*(*sensRsp)[i];
560 Teuchos::RCP<const std::vector<Real> > vp =
562 Teuchos::RCP<const std::vector<Real> > Sp =
567 Real h = std::max(1.0,S.
norm()/v.
norm())*tol;
571 Real Is_scale = pow( 10,
int( log10( (*Sp)[0] ) ) );
572 Real Rs_scale = pow( 10,
int( log10( (*Sp)[1] ) ) );
574 Real h1 = Is_scale*h;
575 Real h2 = Rs_scale*h;
578 Teuchos::RCP<Vector<Real> > g = S.
clone();
582 Teuchos::RCP<std::vector<Real> > Snewp = Teuchos::rcp(
new std::vector<Real> (2, 0.0) );
584 (*Snewp)[0] = (*Sp)[0] + h1*(*vp)[0];
585 (*Snewp)[1] = (*Sp)[1] + h2*(*vp)[1];
593 hv.
scale(1.0/std::sqrt(h1*h1+h2*h2));
596 Teuchos::RCP<std::vector<Real> > hvp =
597 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(hv)).getVector());
598 Teuchos::RCP<const std::vector<Real> > vp =
600 Teuchos::RCP<const std::vector<Real> > Sp =
606 Teuchos::RCP<std::vector<Real> > Ip =
607 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(I)).getVector());
613 StdVector<Real> lambda( Teuchos::rcp(
new std::vector<Real>(n,0.0) ) );
614 Teuchos::RCP<std::vector<Real> > lambdap =
615 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(lambda)).getVector());
621 Teuchos::RCP<std::vector<Real> > wp =
622 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(w)).getVector());
625 for(
int i=0;i<n;i++){
626 (*wp)[i] = ( (*vp)[0] *
diodeIs( (*Ip)[i],(*
Vsrc_)[i],(*Sp)[0],(*Sp)[1] ) + (*vp)[1] *
diodeRs( (*Ip)[i],(*
Vsrc_)[i],(*Sp)[0],(*Sp)[1] ) ) /
diodeI((*Ip)[i],(*Vsrc_)[i],(*Sp)[0],(*Sp)[1]);
630 Teuchos::RCP<std::vector<Real> > pp =
631 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(p)).getVector());
634 for(
int j=0;j<n;j++){
635 (*pp)[j] = ( (*wp)[j] + (*lambdap)[j] *
diodeII( (*Ip)[j],(*
Vsrc_)[j],(*Sp)[0],(*Sp)[1] ) * (*wp)[j] - (*lambdap)[j] *
diodeIIs( (*Ip)[j],(*
Vsrc_)[j],(*Sp)[0],(*Sp)[1] ) * (*vp)[0] - (*lambdap)[j] *
diodeIRs( (*Ip)[j],(*
Vsrc_)[j],(*Sp)[0],(*Sp)[1] ) * (*vp)[1] ) /
diodeI( (*Ip)[j],(*Vsrc_)[j],(*Sp)[0],(*Sp)[1] );
639 (*hvp)[0] = 0.0;(*hvp)[1] = 0.0;
640 for(
int k=0;k<n;k++){
641 (*hvp)[0] +=
diodeIs( (*Ip)[k],(*
Vsrc_)[k],(*Sp)[0],(*Sp)[1] ) * (*pp)[k] - (*lambdap)[k] * (*wp)[k] *
diodeIIs( (*Ip)[k],(*
Vsrc_)[k],(*Sp)[0],(*Sp)[1] ) + (*lambdap)[k] * (*vp)[0] *
diodeIsIs( (*Ip)[k],(*
Vsrc_)[k],(*Sp)[0],(*Sp)[1] ) + (*lambdap)[k] * (*vp)[1] *
diodeIsRs( (*Ip)[k],(*
Vsrc_)[k],(*Sp)[0],(*Sp)[1] );
642 (*hvp)[1] +=
diodeRs( (*Ip)[k],(*
Vsrc_)[k],(*Sp)[0],(*Sp)[1] ) * (*pp)[k] - (*lambdap)[k] * (*wp)[k] *
diodeIRs( (*Ip)[k],(*
Vsrc_)[k],(*Sp)[0],(*Sp)[1] ) + (*lambdap)[k] * (*vp)[0] *
diodeIsRs( (*Ip)[k],(*
Vsrc_)[k],(*Sp)[0],(*Sp)[1] ) + (*lambdap)[k] * (*vp)[1] *
diodeRsRs( (*Ip)[k],(*
Vsrc_)[k],(*Sp)[0],(*Sp)[1] );
647 Teuchos::RCP<std::vector<Real> > hvp =
648 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(hv)).getVector());
649 Teuchos::RCP<const std::vector<Real> > vp =
651 Teuchos::RCP<const std::vector<Real> > Sp =
657 Teuchos::RCP<std::vector<Real> > Ip =
658 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(I)).getVector());
664 StdVector<Real> sensIs( Teuchos::rcp(
new std::vector<Real>(n,0.0) ) );
665 StdVector<Real> sensRs( Teuchos::rcp(
new std::vector<Real>(n,0.0) ) );
669 Teuchos::RCP<std::vector<Real> > sensIsp = Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(sensIs)).getVector());
670 Teuchos::RCP<std::vector<Real> > sensRsp = Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(sensRs)).getVector());
673 Real H11 = 0.0; Real H12 = 0.0; Real H22 = 0.0;
674 for(
int k=0;k<n;k++){
675 H11 += (*sensIsp)[k]*(*sensIsp)[k];
676 H12 += (*sensIsp)[k]*(*sensRsp)[k];
677 H22 += (*sensRsp)[k]*(*sensRsp)[k];
681 (*hvp)[0] = H11*(*vp)[0] + H12*(*vp)[1];
682 (*hvp)[1] = H12*(*vp)[0] + H22*(*vp)[1];
700 void generate_plot(Real Is_lo, Real Is_up, Real Is_step, Real Rs_lo, Real Rs_up, Real Rs_step){
701 Teuchos::RCP<std::vector<double> > S_rcp = Teuchos::rcp(
new std::vector<double>(2,0.0) );
703 std::ofstream output (
"Objective.dat");
709 int n = (Is_up-Is_lo)/Is_step + 1;
710 int m = (Rs_up-Rs_lo)/Rs_step + 1;
711 for(
int i=0;i<n;i++){
712 Is = Is_lo + i*Is_step;
713 for(
int j=0;j<m;j++){
714 Rs = Rs_lo + j*Rs_step;
717 val = this->
value(S,tol);
718 if(output.is_open()){
719 output << std::scientific << Is <<
" " << Rs <<
" " << val <<
"\n";
748 x_lo_.push_back(lo_Is);
749 x_lo_.push_back(lo_Rs);
751 x_up_.push_back(up_Is);
752 x_up_.push_back(up_Rs);
758 Teuchos::RCP<std::vector<Real> > ex =
759 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(x)).getVector());
760 (*ex)[0] = std::max(
x_lo_[0],std::min(
x_up_[0],(*ex)[0]));
761 (*ex)[1] = std::max(
x_lo_[1],std::min(
x_up_[1],(*ex)[1]));
764 Teuchos::RCP<const std::vector<Real> > ex =
766 return ((*ex)[0] >= this->
x_lo_[0] && (*ex)[1] >= this->
x_lo_[1] &&
767 (*ex)[0] <= this->
x_up_[0] && (*ex)[1] <= this->
x_up_[1]);
770 Teuchos::RCP<const std::vector<Real> > ex =
772 Teuchos::RCP<std::vector<Real> > ev =
773 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(v)).getVector());
776 for (
int i = 0; i < 2; i++ ) {
777 if ( ((*ex)[i] <= this->
x_lo_[i]+epsn) ) {
783 Teuchos::RCP<const std::vector<Real> > ex =
785 Teuchos::RCP<std::vector<Real> > ev =
786 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(v)).getVector());
789 for (
int i = 0; i < 2; i++ ) {
790 if ( ((*ex)[i] >= this->
x_up_[i]-epsn) ) {
796 Teuchos::RCP<const std::vector<Real> > ex =
798 Teuchos::RCP<std::vector<Real> > ev =
799 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(v)).getVector());
802 for (
int i = 0; i < 2; i++ ) {
803 if ( ((*ex)[i] <= this->
x_lo_[i]+epsn) ||
804 ((*ex)[i] >= this->
x_up_[i]-epsn) ) {
810 Teuchos::RCP<const std::vector<Real> > ex =
812 Teuchos::RCP<const std::vector<Real> > eg =
814 Teuchos::RCP<std::vector<Real> > ev =
815 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(v)).getVector());
818 for (
int i = 0; i < 2; i++ ) {
819 if ( ((*ex)[i] <= this->
x_lo_[i]+epsn && (*eg)[i] > 0.0) ) {
825 Teuchos::RCP<const std::vector<Real> > ex =
827 Teuchos::RCP<const std::vector<Real> > eg =
829 Teuchos::RCP<std::vector<Real> > ev =
830 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(v)).getVector());
833 for (
int i = 0; i < 2; i++ ) {
834 if ( ((*ex)[i] >= this->
x_up_[i]-epsn && (*eg)[i] < 0.0) ) {
840 Teuchos::RCP<const std::vector<Real> > ex =
842 Teuchos::RCP<const std::vector<Real> > eg =
844 Teuchos::RCP<std::vector<Real> > ev =
845 Teuchos::rcp_const_cast<std::vector<Real> >((Teuchos::dyn_cast<
StdVector<Real> >(v)).getVector());
848 for (
int i = 0; i < 2; i++ ) {
849 if ( ((*ex)[i] <= this->
x_lo_[i]+epsn && (*eg)[i] > 0.0) ||
850 ((*ex)[i] >= this->
x_up_[i]-epsn && (*eg)[i] < 0.0) ) {
857 Teuchos::RCP<std::vector<Real> > us = Teuchos::rcp(
new std::vector<Real>(2,0.0) );
858 us->assign(this->
x_up_.begin(),this->
x_up_.end());
864 Teuchos::RCP<std::vector<Real> > ls = Teuchos::rcp(
new std::vector<Real>(2,0.0) );
865 ls->assign(this->
x_lo_.begin(),this->
x_lo_.end());
Provides the interface to evaluate objective functions.
bool use_adjoint_
If true, use adjoint gradient computation, else compute gradient using sensitivities.
Real diodeIsRs(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Second derivative of diode equation wrt Is and Rs.
Real noise_
Percentage of noise to add to measurements; if 0.0 - no noise.
Real diodeIsIs(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Second derivative of diode equation wrt Is^2.
void pruneUpperActive(Vector< Real > &v, const Vector< Real > &x, Real eps)
Set variables to zero if they correspond to the upper -active set.
virtual void scale(const Real alpha)=0
Compute where .
Real min_diff_
Half of the minimum distance between upper and lower bounds.
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
Objective_DiodeCircuit(Real Vth, std::ifstream &input_file, bool lambertw, Real noise, bool use_adjoint, int use_hessvec)
A constructor using data from given file.
Real scale_
Scaling for the epsilon margin.
BoundConstraint_DiodeCircuit(Real scale, Real lo_Is, Real up_Is, Real lo_Rs, Real up_Rs)
Real lambertWCurrent(Real Is, Real Rs, Real Vsrc)
Find currents using Lambert-W function.
int use_hessvec_
0 - use FD(with scaling), 1 - use exact implementation (with second order derivatives), 2 - use Gauss-Newton approximation (first order derivatives only)
void pruneActive(Vector< Real > &v, const Vector< Real > &x, Real eps)
Set variables to zero if they correspond to the -active set.
virtual void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
Apply Hessian approximation to vector.
void solve_sensitivity_Is(Vector< Real > &sens, const Vector< Real > &I, const Vector< Real > &S)
Solve the sensitivity equation wrt Is.
void solve_adjoint(Vector< Real > &lambda, const Vector< Real > &I, const Vector< Real > &S)
Solve the adjoint equation.
void gradient(Vector< Real > &g, const Vector< Real > &S, Real &tol)
Compute the gradient of the reduced objective function either using adjoint or using sensitivities...
Teuchos::RCP< std::vector< Real > > Imeas_
Vector of measured currents in DC analysis (data)
Real Vth_
Thermal voltage (constant)
virtual Teuchos::RCP< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
void setVectorToUpperBound(ROL::Vector< Real > &u)
Set the input vector to the upper bound.
void generate_plot(Real Is_lo, Real Is_up, Real Is_step, Real Rs_lo, Real Rs_up, Real Rs_step)
Generate data to plot objective function.
virtual void zero()
Set to zero vector.
Defines the linear algebra or vector space interface.
Bound constraints on optimization parameters.
Teuchos::RCP< std::vector< Real > > Vsrc_
Vector of source voltages in DC analysis (input)
void pruneLowerActive(Vector< Real > &v, const Vector< Real > &x, Real eps)
Set variables to zero if they correspond to the lower -active set.
The diode circuit problem.
Real value(const Vector< Real > &S, Real &tol)
Evaluate objective function.
void pruneLowerActive(Vector< Real > &v, const Vector< Real > &g, const Vector< Real > &x, Real eps)
Set variables to zero if they correspond to the lower -binding set.
Real diodeI(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Derivative of diode equation wrt I.
Provides the std::vector implementation of the ROL::Vector interface.
std::vector< Real > x_lo_
Vector of lower bounds.
void pruneActive(Vector< Real > &v, const Vector< Real > &g, const Vector< Real > &x, Real eps)
Set variables to zero if they correspond to the -binding set.
bool isFeasible(const Vector< Real > &x)
Check if the vector, v, is feasible.
void lambertw(Real x, Real &w, int &ierr, Real &xi)
Lambert-W function for diodes.
Real Newton(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Newton's method with line search.
void hessVec(Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &S, Real &tol)
Compute the Hessian-vector product of the reduced objective function.
Provides the interface to apply upper and lower bound constraints.
void solve_circuit(Vector< Real > &I, const Vector< Real > &S)
Solve circuit given optimization parameters Is and Rs.
void project(Vector< Real > &x)
Project optimization variables onto the bounds.
std::vector< Real > x_up_
Vector of upper bounds.
bool lambertw_
If true, use Lambert-W function to solve circuit, else use Newton's method.
virtual void set(const Vector &x)
Set where .
virtual Real norm() const =0
Returns where .
Real diode(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Diode equation.
Real diodeIs(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Derivative of diode equation wrt Is.
Real diodeRsRs(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Second derivative of diode equation wrt Rs^2.
Objective_DiodeCircuit(Real Vth, Real Vsrc_min, Real Vsrc_max, Real Vsrc_step, Real true_Is, Real true_Rs, bool lambertw, Real noise, bool use_adjoint, int use_hessvec)
A constructor generating data.
void solve_sensitivity_Rs(Vector< Real > &sens, const Vector< Real > &I, const Vector< Real > &S)
Solve the sensitivity equation wrt Rs.
void setVectorToLowerBound(ROL::Vector< Real > &l)
Set the input vector to the lower bound.
void pruneUpperActive(Vector< Real > &v, const Vector< Real > &g, const Vector< Real > &x, Real eps)
Set variables to zero if they correspond to the upper -binding set.
Real diodeIIs(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Second derivative of diode equation wrt I and Is.
Real diodeIRs(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Second derivative of diode equation wrt I and Rs.
Real diodeII(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Second derivative of diode equation wrt I^2.
static const double ROL_EPSILON
Platform-dependent machine epsilon.
Real diodeRs(const Real I, const Real Vsrc, const Real Is, const Real Rs)
Derivative of diode equation wrt Rs.
void set_method(bool lambertw)
Change the method for solving the circuit if needed.