ROL
Public Member Functions | Private Attributes | List of all members
ROL::Objective_DiodeCircuit< Real > Class Template Reference

The diode circuit problem. More...

#include <ROL_DiodeCircuit.hpp>

+ Inheritance diagram for ROL::Objective_DiodeCircuit< Real >:

Public Member Functions

 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. More...
 
 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. More...
 
void set_method (bool lambertw)
 Change the method for solving the circuit if needed.
 
Real diode (const Real I, const Real Vsrc, const Real Is, const Real Rs)
 Diode equation. More...
 
Real diodeI (const Real I, const Real Vsrc, const Real Is, const Real Rs)
 Derivative of diode equation wrt I.
 
Real diodeIs (const Real I, const Real Vsrc, const Real Is, const Real Rs)
 Derivative of diode equation wrt Is.
 
Real diodeRs (const Real I, const Real Vsrc, const Real Is, const Real Rs)
 Derivative of diode equation wrt Rs.
 
Real diodeII (const Real I, const Real Vsrc, const Real Is, const Real Rs)
 Second derivative of diode equation wrt I^2.
 
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 diodeIsIs (const Real I, const Real Vsrc, const Real Is, const Real Rs)
 Second derivative of diode equation wrt Is^2.
 
Real diodeIsRs (const Real I, const Real Vsrc, const Real Is, const Real Rs)
 Second derivative of diode equation wrt Is and Rs.
 
Real diodeRsRs (const Real I, const Real Vsrc, const Real Is, const Real Rs)
 Second derivative of diode equation wrt Rs^2.
 
Real Newton (const Real I, const Real Vsrc, const Real Is, const Real Rs)
 Newton's method with line search. More...
 
void lambertw (Real x, Real &w, int &ierr, Real &xi)
 Lambert-W function for diodes. More...
 
Real lambertWCurrent (Real Is, Real Rs, Real Vsrc)
 Find currents using Lambert-W function. More...
 
void solve_circuit (Vector< Real > &I, const Vector< Real > &S)
 Solve circuit given optimization parameters Is and Rs.
 
Real value (const Vector< Real > &S, Real &tol)
 Evaluate objective function. More...
 
void solve_adjoint (Vector< Real > &lambda, const Vector< Real > &I, const Vector< Real > &S)
 Solve the adjoint equation. More...
 
void solve_sensitivity_Is (Vector< Real > &sens, const Vector< Real > &I, const Vector< Real > &S)
 Solve the sensitivity equation wrt Is. More...
 
void solve_sensitivity_Rs (Vector< Real > &sens, const Vector< Real > &I, const Vector< Real > &S)
 Solve the sensitivity equation wrt Rs. More...
 
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.
 
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. More...
 
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. More...
 
- Public Member Functions inherited from ROL::Objective< Real >
virtual void update (const Vector< Real > &x, bool flag=true, int iter=-1)
 Update objective function. More...
 
virtual Real dirDeriv (const Vector< Real > &x, const Vector< Real > &d, Real &tol)
 Compute directional derivative. More...
 
virtual void invHessVec (Vector< Real > &hv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
 Apply inverse Hessian approximation to vector. More...
 
virtual void precond (Vector< Real > &Pv, const Vector< Real > &v, const Vector< Real > &x, Real &tol)
 Apply preconditioner to vector. More...
 
virtual std::vector
< std::vector< Real > > 
checkGradient (const Vector< Real > &x, const Vector< Real > &d, const bool printToScreen=true, const int numSteps=ROL_NUM_CHECKDERIV_STEPS)
 Finite-difference gradient check. More...
 
virtual std::vector
< std::vector< Real > > 
checkHessVec (const Vector< Real > &x, const Vector< Real > &v, const bool printToScreen=true, const int numSteps=ROL_NUM_CHECKDERIV_STEPS)
 Finite-difference Hessian-applied-to-vector check. More...
 
virtual std::vector< Real > checkHessSym (const Vector< Real > &x, const Vector< Real > &v, const Vector< Real > &w, const bool printToScreen=true)
 Hessian symmetry check. More...
 

Private Attributes

Real Vth_
 Thermal voltage (constant)
 
Teuchos::RCP< std::vector< Real > > Imeas_
 Vector of measured currents in DC analysis (data)
 
Teuchos::RCP< std::vector< Real > > Vsrc_
 Vector of source voltages in DC analysis (input)
 
bool lambertw_
 If true, use Lambert-W function to solve circuit, else use Newton's method.
 
Real noise_
 Percentage of noise to add to measurements; if 0.0 - no noise.
 
bool use_adjoint_
 If true, use adjoint gradient computation, else compute gradient using sensitivities.
 
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)
 

Detailed Description

template<class Real>
class ROL::Objective_DiodeCircuit< Real >

The diode circuit problem.

The diode circuit problem:

\begin{eqnarray*} \min_{I_S,R_S} \,\, \frac{1}{2}\sum\limits_{n=1}^N (I_n-I_n^{meas})^2 \\ \text{s.t.}\;\;\begin{cases}c(I_S,R_S,I_1,V^{src}_1)=0\\ \dots \\c(I_S,R_S,I_N,V^{src}_N)=0\end{cases} \end{eqnarray*}

where

\[c(I_S,R_S,I_n,V^{src}_n)=I_n - I_S\left(\exp\left(\frac{-I_n R_S+V^{src}_n}{V_{th}}\right)-1\right)\]

.

Definition at line 32 of file ROL_DiodeCircuit.hpp.

Constructor & Destructor Documentation

template<class Real >
ROL::Objective_DiodeCircuit< Real >::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 
)
inline

A constructor generating data.

Given thermal voltage, minimum and maximum values of source voltages and a step size, values of Is and Rs generates vector of source voltages and solves nonlinear diode equation to populate the vector of measured currents, which is later used as data. If noise is nonzero, adds random perturbation to data on the order of the magnitude of the components. Sets the flag to use Lambert-W function or Newton's method to solve circuit. Sets the flags to use adjoint gradient computation and one of three Hessian-vector implementations.


Definition at line 58 of file ROL_DiodeCircuit.hpp.

References ROL::Objective_DiodeCircuit< Real >::Imeas_, ROL::Objective_DiodeCircuit< Real >::lambertw(), ROL::Objective_DiodeCircuit< Real >::lambertw_, ROL::Objective_DiodeCircuit< Real >::lambertWCurrent(), ROL::Objective_DiodeCircuit< Real >::Newton(), ROL::Objective_DiodeCircuit< Real >::use_adjoint_, ROL::Objective_DiodeCircuit< Real >::use_hessvec_, ROL::Objective_DiodeCircuit< Real >::Vsrc_, and ROL::Objective_DiodeCircuit< Real >::Vth_.

template<class Real >
ROL::Objective_DiodeCircuit< Real >::Objective_DiodeCircuit ( Real  Vth,
std::ifstream &  input_file,
bool  lambertw,
Real  noise,
bool  use_adjoint,
int  use_hessvec 
)
inline

A constructor using data from given file.

Given thermal voltage and a file with two columns - one for source voltages, another for corresponding currents - populates vectors of source voltages and measured currents. If noise is nonzero, adds random perturbation to data on the order of the magnitude of the components. Sets the flag to use Lambert-W function or Newton's method to solve circuit. Sets the flags to use adjoint gradient computation and one of three Hessian-vector implementations.


Definition at line 110 of file ROL_DiodeCircuit.hpp.

References ROL::Objective_DiodeCircuit< Real >::Imeas_, ROL::Objective_DiodeCircuit< Real >::lambertw(), ROL::Objective_DiodeCircuit< Real >::lambertw_, ROL::Objective_DiodeCircuit< Real >::use_adjoint_, ROL::Objective_DiodeCircuit< Real >::use_hessvec_, ROL::Objective_DiodeCircuit< Real >::Vsrc_, and ROL::Objective_DiodeCircuit< Real >::Vth_.

Member Function Documentation

template<class Real >
Real ROL::Objective_DiodeCircuit< Real >::diode ( const Real  I,
const Real  Vsrc,
const Real  Is,
const Real  Rs 
)
inline

Diode equation.

Diode equation formula: \( I-I_S\left(\exp\left(\frac{V_{src}-IR_S}{V_{th}}\right)-1\right) \).


Definition at line 148 of file ROL_DiodeCircuit.hpp.

References ROL::Objective_DiodeCircuit< Real >::Vth_.

Referenced by ROL::Objective_DiodeCircuit< Real >::Newton().

template<class Real >
Real ROL::Objective_DiodeCircuit< Real >::Newton ( const Real  I,
const Real  Vsrc,
const Real  Is,
const Real  Rs 
)
inline

Newton's method with line search.

Solves the diode equation for the current using Newton's method.


Definition at line 204 of file ROL_DiodeCircuit.hpp.

References ROL::Objective_DiodeCircuit< Real >::diode(), and ROL::Objective_DiodeCircuit< Real >::diodeI().

Referenced by ROL::Objective_DiodeCircuit< Real >::Objective_DiodeCircuit(), and ROL::Objective_DiodeCircuit< Real >::solve_circuit().

template<class Real >
void ROL::Objective_DiodeCircuit< Real >::lambertw ( Real  x,
Real &  w,
int &  ierr,
Real &  xi 
)
inline

Lambert-W function for diodes.

Function : DeviceSupport::lambertw Purpose : provides a lambert-w function for diodes and BJT's. Special Notes :

Purpose. Evaluate principal branch of Lambert W function at x.

w = w(x) is the value of Lambert's function. ierr = 0 indicates a safe return. ierr = 1 if x is not in the domain. ierr = 2 if the computer arithmetic contains a bug. xi may be disregarded (it is the error).

Prototype: void lambertw( double, double, int, double);

Reference: T.C. Banwell Bipolar transistor circuit analysis using the Lambert W-function, IEEE Transactions on Circuits and Systems I: Fundamental Theory and Applications

vol. 47, pp. 1621-1633, Nov. 2000.

Scope : public Creator : David Day, SNL Creation Date : 04/16/02


Definition at line 277 of file ROL_DiodeCircuit.hpp.

Referenced by ROL::Objective_DiodeCircuit< Real >::lambertWCurrent(), ROL::Objective_DiodeCircuit< Real >::Objective_DiodeCircuit(), and ROL::Objective_DiodeCircuit< Real >::set_method().

template<class Real >
Real ROL::Objective_DiodeCircuit< Real >::lambertWCurrent ( Real  Is,
Real  Rs,
Real  Vsrc 
)
inline

Find currents using Lambert-W function.

Reference: T.C. Banwell Bipolar transistor circuit analysis using the Lambert W-function, IEEE Transactions on Circuits and Systems I: Fundamental Theory and Applications vol. 47, pp. 1621-1633, Nov. 2000.


Definition at line 350 of file ROL_DiodeCircuit.hpp.

References ROL::Objective_DiodeCircuit< Real >::lambertw(), and ROL::Objective_DiodeCircuit< Real >::Vth_.

Referenced by ROL::Objective_DiodeCircuit< Real >::Objective_DiodeCircuit(), and ROL::Objective_DiodeCircuit< Real >::solve_circuit().

template<class Real >
Real ROL::Objective_DiodeCircuit< Real >::value ( const Vector< Real > &  S,
Real &  tol 
)
inlinevirtual

Evaluate objective function.

\(\frac{1}{2}\sum\limits_{i=1}^{N}(I_i-I^{meas}_i)^2\)


Implements ROL::Objective< Real >.

Definition at line 400 of file ROL_DiodeCircuit.hpp.

References ROL::Objective_DiodeCircuit< Real >::Imeas_, and ROL::Objective_DiodeCircuit< Real >::solve_circuit().

Referenced by ROL::Objective_DiodeCircuit< Real >::generate_plot().

template<class Real >
void ROL::Objective_DiodeCircuit< Real >::solve_adjoint ( Vector< Real > &  lambda,
const Vector< Real > &  I,
const Vector< Real > &  S 
)
inline

Solve the adjoint equation.

\(\lambda_i = \frac{(I^{meas}_i-I_i)}{\frac{\partial c}{\partial I}(I_i,V^{src}_i,I_S,R_S)}\)


Definition at line 425 of file ROL_DiodeCircuit.hpp.

References ROL::Objective_DiodeCircuit< Real >::diodeI().

Referenced by ROL::Objective_DiodeCircuit< Real >::gradient(), and ROL::Objective_DiodeCircuit< Real >::hessVec().

template<class Real >
void ROL::Objective_DiodeCircuit< Real >::solve_sensitivity_Is ( Vector< Real > &  sens,
const Vector< Real > &  I,
const Vector< Real > &  S 
)
inline

Solve the sensitivity equation wrt Is.

Computes sensitivity

\[\frac{\partial I}{\partial Is}\]


Definition at line 448 of file ROL_DiodeCircuit.hpp.

References ROL::Objective_DiodeCircuit< Real >::diodeI(), ROL::Objective_DiodeCircuit< Real >::diodeIs(), and ROL::Objective_DiodeCircuit< Real >::Vsrc_.

Referenced by ROL::Objective_DiodeCircuit< Real >::gradient(), and ROL::Objective_DiodeCircuit< Real >::hessVec().

template<class Real >
void ROL::Objective_DiodeCircuit< Real >::solve_sensitivity_Rs ( Vector< Real > &  sens,
const Vector< Real > &  I,
const Vector< Real > &  S 
)
inline

Solve the sensitivity equation wrt Rs.

Computes sensitivity

\[\frac{\partial I}{\partial Rs}\]


Definition at line 470 of file ROL_DiodeCircuit.hpp.

References ROL::Objective_DiodeCircuit< Real >::diodeI(), ROL::Objective_DiodeCircuit< Real >::diodeRs(), and ROL::Objective_DiodeCircuit< Real >::Vsrc_.

Referenced by ROL::Objective_DiodeCircuit< Real >::gradient(), and ROL::Objective_DiodeCircuit< Real >::hessVec().

template<class Real >
void ROL::Objective_DiodeCircuit< Real >::hessVec ( Vector< Real > &  hv,
const Vector< Real > &  v,
const Vector< Real > &  S,
Real &  tol 
)
inlinevirtual
template<class Real >
void ROL::Objective_DiodeCircuit< Real >::generate_plot ( Real  Is_lo,
Real  Is_up,
Real  Is_step,
Real  Rs_lo,
Real  Rs_up,
Real  Rs_step 
)
inline

Generate data to plot objective function.

Generates a file with three columns - Is value, Rs value, objective value. To plot with gnuplot type: gnuplot; set dgrid3d 100,100; set hidden3d; splot "Objective.dat" u 1:2:3 with lines;


Definition at line 699 of file ROL_DiodeCircuit.hpp.

References ROL::Objective_DiodeCircuit< Real >::value().


The documentation for this class was generated from the following file: