ROL
|
The diode circuit problem. More...
#include <ROL_DiodeCircuit.hpp>
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... | |
![]() | |
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) | |
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.
|
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_.
|
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_.
|
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().
|
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().
|
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().
|
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().
|
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().
|
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().
|
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().
|
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().
|
inlinevirtual |
Compute the Hessian-vector product of the reduced objective function.
Hessian-times-vector computation.
Reimplemented from ROL::Objective< Real >.
Definition at line 555 of file ROL_DiodeCircuit.hpp.
References ROL::Vector< Real >::axpy(), ROL::Vector< Real >::clone(), ROL::Objective_DiodeCircuit< Real >::diodeI(), ROL::Objective_DiodeCircuit< Real >::diodeII(), ROL::Objective_DiodeCircuit< Real >::diodeIIs(), ROL::Objective_DiodeCircuit< Real >::diodeIRs(), ROL::Objective_DiodeCircuit< Real >::diodeIs(), ROL::Objective_DiodeCircuit< Real >::diodeIsIs(), ROL::Objective_DiodeCircuit< Real >::diodeIsRs(), ROL::Objective_DiodeCircuit< Real >::diodeRs(), ROL::Objective_DiodeCircuit< Real >::diodeRsRs(), ROL::Objective_DiodeCircuit< Real >::gradient(), ROL::Objective< Real >::hessVec(), ROL::Objective_DiodeCircuit< Real >::Imeas_, ROL::Vector< Real >::norm(), ROL::Vector< Real >::scale(), ROL::Objective_DiodeCircuit< Real >::solve_adjoint(), ROL::Objective_DiodeCircuit< Real >::solve_circuit(), ROL::Objective_DiodeCircuit< Real >::solve_sensitivity_Is(), ROL::Objective_DiodeCircuit< Real >::solve_sensitivity_Rs(), ROL::Objective_DiodeCircuit< Real >::use_hessvec_, ROL::Objective_DiodeCircuit< Real >::Vsrc_, and ROL::Vector< Real >::zero().
|
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().