twoproj.c

来自「ngspice又一个电子CAD仿真软件代码.功能更全」· C语言 代码 · 共 680 行 · 第 1/2 页

/**********
Copyright 1991 Regents of the University of California.  All rights reserved.
Author:	1987 Kartikeya Mayaram, U. C. Berkeley CAD Group
Author:	1991 David A. Gates, U. C. Berkeley CAD Group
**********/

/*
 * Functions for projecting the next solution point for use with the modified 
 * two-level Newton scheme
 */

#include "ngspice.h"
#include "numglobs.h"
#include "numenum.h"
#include "twomesh.h"
#include "twodev.h"
#include "bool.h"
#include "spmatrix.h"
#include "twoddefs.h"
#include "twodext.h"
#include "cidersupt.h"


/* Forward Declarations */


void NUMD2project(TWOdevice *pDevice, double delV)
{
  TWOnode *pNode;
  TWOelem *pElem;
  int index, eIndex, numContactNodes;
  TWOcontact *pContact = pDevice->pLastContact;
  double *incVpn, *solution = pDevice->dcSolution;
  double delPsi, delN, delP, newN, newP;
  
  delV = -delV / VNorm;
  /* update the boundary condition for the last contact */
  numContactNodes = pContact->numNodes;
  for ( index = 0; index < numContactNodes; index++ ) {
    pNode = pContact->pNodes[ index ];
    pNode->psi += delV;
  }
  
  /* 
   * store the new rhs for computing the incremental quantities
   * with the second to last node. solve the system of equations
   */
  
  if ( ABS(delV) < MIN_DELV ) {
    TWOstoreInitialGuess( pDevice );
    return;
  }
  incVpn = pDevice->dcDeltaSolution;
  storeNewRhs( pDevice, pContact );
  spSolve( pDevice->matrix, pDevice->rhs, incVpn, NIL(spREAL), NIL(spREAL) );
  
  for ( eIndex = 1; eIndex <= pDevice->numElems; eIndex++ ) {
    pElem = pDevice->elements[ eIndex ];
    for ( index = 0; index <= 3; index++ ) {
      if ( pElem->evalNodes[ index ] ) {
	pNode = pElem->pNodes[ index ];
	if ( pNode->nodeType != CONTACT ) {
	  delPsi = incVpn[ pNode->psiEqn ] * delV;
	  solution[ pNode->psiEqn ] = pNode->psi + delPsi;
	  if ( pElem->elemType == SEMICON
	      && (!OneCarrier || (OneCarrier == N_TYPE)) ) {
	    delN = incVpn[ pNode->nEqn ] * delV;
	    newN = pNode->nConc + delN;
	    if ( newN <= 0.0 ) {
	      solution[ pNode->nEqn ] = guessNewConc( pNode->nConc, delN );
	    }
	    else {
	      solution[ pNode->nEqn ] = newN;
	    }
	  }
	  if ( pElem->elemType == SEMICON
	      && (!OneCarrier || (OneCarrier == P_TYPE)) ) {
	    delP = incVpn[ pNode->pEqn ] * delV;
	    newP = pNode->pConc + delP;
	    if ( newP <= 0.0 ) {
	      solution[ pNode->pEqn ] = guessNewConc( pNode->pConc, delP );
	    }
	    else {
	      solution[ pNode->pEqn ] = newP;
	    }
	  }
	}
      }
    }
  }
}


void NBJT2project(TWOdevice *pDevice, double delVce, double delVbe)
{
  TWOnode *pNode;
  TWOelem *pElem;
  int index, eIndex, numContactNodes;
  TWOcontact *pColContact = pDevice->pFirstContact;
  TWOcontact *pBaseContact = pDevice->pFirstContact->next;
  double *incVce, *incVbe, *solution = pDevice->dcSolution;
  double delPsi, delN, delP, newN, newP;
  double nConc, pConc;
  
  /* Normalize the voltages for calculations. */
  if ( delVce != 0.0 ) {
    delVce = delVce / VNorm;
    numContactNodes = pColContact->numNodes;
    for ( index = 0; index < numContactNodes; index++ ) {
      pNode = pColContact->pNodes[ index ];
      pNode->psi += delVce;
    }
  }
  if ( delVbe != 0.0 ) {
    delVbe = delVbe / VNorm;
    numContactNodes = pBaseContact->numNodes;
    for ( index = 0; index < numContactNodes; index++ ) {
      pNode = pBaseContact->pNodes[ index ];
      pNode->psi += delVbe;
    }
  }
  
  /* 
   * store the new rhs for computing the incremental quantities
   * incVce (dcDeltaSolution) and incVbe (copiedSolution) are used to
   * store the incremental quantities associated with Vce and Vbe
   */
  
  /* set incVce = dcDeltaSolution; incVbe = copiedSolution */
  
  if ( ABS( delVce ) > MIN_DELV ) {
    incVce = pDevice->dcDeltaSolution;
    storeNewRhs( pDevice, pColContact );
    spSolve( pDevice->matrix, pDevice->rhs, incVce, NIL(spREAL), NIL(spREAL));
    
    for ( eIndex = 1; eIndex <= pDevice->numElems; eIndex++ ) {
      pElem = pDevice->elements[ eIndex ];
      for ( index = 0; index <= 3; index++ ) {
	if ( pElem->evalNodes[ index ] ) {
	  pNode = pElem->pNodes[ index ];
	  if ( pNode->nodeType != CONTACT ) {
	    delPsi = incVce[ pNode->psiEqn ] * delVce;
	    solution[ pNode->psiEqn ] = pNode->psi + delPsi;
	    if ( pElem->elemType == SEMICON
		&& (!OneCarrier || (OneCarrier == N_TYPE)) ) {
	      delN = incVce[ pNode->nEqn ] * delVce;
	      newN = pNode->nConc + delN;
	      if ( newN <= 0.0 ) {
		solution[ pNode->nEqn ] = guessNewConc( pNode->nConc, delN );
	      }
	      else {
		solution[ pNode->nEqn ] = newN;
	      }
	    }
	    if ( pElem->elemType == SEMICON
		&& (!OneCarrier | (OneCarrier == P_TYPE)) ) {
	      delP = incVce[ pNode->pEqn ] * delVce;
	      newP = pNode->pConc + delP;
	      if ( newP <= 0.0 ) {
		solution[ pNode->pEqn ] = guessNewConc( pNode->pConc, delP );
	      }
	      else {
		solution[ pNode->pEqn ] = newP;
	      }
	    }
	  }
	}
      }
    }
  }
  else {
    TWOstoreInitialGuess( pDevice );
  }
  
  if ( ABS( delVbe ) > MIN_DELV ) {
    incVbe = pDevice->copiedSolution;
    storeNewRhs( pDevice, pBaseContact );
    spSolve( pDevice->matrix, pDevice->rhs, incVbe, NIL(spREAL), NIL(spREAL));
    
    for ( eIndex = 1; eIndex <= pDevice->numElems; eIndex++ ) {
      pElem = pDevice->elements[ eIndex ];
      for ( index = 0; index <= 3; index++ ) {
	if ( pElem->evalNodes[ index ] ) {
	  pNode = pElem->pNodes[ index ];
	  if ( pNode->nodeType != CONTACT ) {
	    delPsi = incVbe[ pNode->psiEqn ] * delVbe;
	    solution[ pNode->psiEqn ] += delPsi;
	    if ( pElem->elemType == SEMICON
		&& (!OneCarrier || (OneCarrier == N_TYPE)) ) {
	      delN = incVbe[ pNode->nEqn ] * delVbe;
	      nConc = solution[ pNode->nEqn ];
	      newN = nConc + delN;
	      if ( newN <= 0.0 ) {
		solution[ pNode->nEqn ] = guessNewConc( nConc, delN );
	      }
	      else {
		solution[ pNode->nEqn ] = newN;
	      }
	    }
	    if ( pElem->elemType == SEMICON
		&& (!OneCarrier || (OneCarrier == P_TYPE)) ) {
	      delP = incVbe[ pNode->pEqn ] * delVbe;
	      pConc = solution[ pNode->pEqn ];
	      newP = pConc + delP;
	      if ( newP <= 0.0 ) {
		solution[ pNode->pEqn ] = guessNewConc( pConc, delP );
	      }
	      else {
		solution[ pNode->pEqn ] = newP;
	      }
	    }
	  }
	}
      }
    }
  }
}


void NUMOSproject(TWOdevice *pDevice, double delVdb, double delVsb, 
                  double delVgb)
{
  TWOnode *pNode;
  TWOelem *pElem;
  int index, eIndex, numContactNodes;
  TWOcontact *pDContact = pDevice->pFirstContact;
  TWOcontact *pGContact = pDevice->pFirstContact->next;
  TWOcontact *pSContact = pDevice->pFirstContact->next->next;
  double *incVdb, *incVsb, *incVgb, *solution = pDevice->dcSolution;
  double delPsi, delN, delP, newN, newP;
  double nConc, pConc;
  
  /* normalize the voltages for calculations */
  if ( delVdb != 0.0 ) {
    delVdb = delVdb / VNorm;
    numContactNodes = pDContact->numNodes;
    for ( index = 0; index < numContactNodes; index++ ) {
      pNode = pDContact->pNodes[ index ];
      pNode->psi += delVdb;
    }
  }
  if ( delVsb != 0.0 ) {
    delVsb = delVsb / VNorm;
    numContactNodes = pSContact->numNodes;
    for ( index = 0; index < numContactNodes; index++ ) {
      pNode = pSContact->pNodes[ index ];
      pNode->psi += delVsb;
    }
  }
  if ( delVgb != 0.0 ) {
    delVgb = delVgb / VNorm;
    numContactNodes = pGContact->numNodes;
    for ( index = 0; index < numContactNodes; index++ ) {
      pNode = pGContact->pNodes[ index ];
      pNode->psi += delVgb;
    }
  }
  
  /* 
   * store the new rhs for computing the incremental quantities
   * incVdb (dcDeltaSolution), incVsb, incVgb 
   */
  
  if ( ABS( delVdb ) > MIN_DELV ) {
    
    incVdb = pDevice->dcDeltaSolution;
    storeNewRhs( pDevice, pDContact );
    spSolve( pDevice->matrix, pDevice->rhs, incVdb, NIL(spREAL), NIL(spREAL));
    
    for ( eIndex = 1; eIndex <= pDevice->numElems; eIndex++ ) {
      pElem = pDevice->elements[ eIndex ];
      for ( index = 0; index <= 3; index++ ) {
	if ( pElem->evalNodes[ index ] ) {
	  pNode = pElem->pNodes[ index ];
	  if ( pNode->nodeType != CONTACT ) {
	    delPsi = incVdb[ pNode->psiEqn ] * delVdb;
	    solution[ pNode->psiEqn ] = pNode->psi + delPsi;
	    if ( pElem->elemType == SEMICON
		&& (!OneCarrier || (OneCarrier == N_TYPE)) ) {
	      delN = incVdb[ pNode->nEqn ] * delVdb;
	      newN = pNode->nConc + delN;
	      if ( newN <= 0.0 ) {
		solution[ pNode->nEqn ] = guessNewConc( pNode->nConc, delN );
	      }
	      else {
		solution[ pNode->nEqn ] = newN;
	      }
	    }
	    if ( pElem->elemType == SEMICON
		&& (!OneCarrier || (OneCarrier == P_TYPE)) ) {
	      delP = incVdb[ pNode->pEqn ] * delVdb;
	      newP = pNode->pConc + delP;
	      if ( newP <= 0.0 ) {
		solution[ pNode->pEqn ] = guessNewConc( pNode->pConc, delP );
	      }
	      else {
		solution[ pNode->pEqn ] = newP;
	      }
	    }
	  }
	}
      }
    }
  }
  else {
    TWOstoreInitialGuess( pDevice );
  }
  
  if ( ABS( delVsb ) > MIN_DELV ) {
    
    incVsb = pDevice->dcDeltaSolution;
    storeNewRhs( pDevice, pSContact );
    spSolve( pDevice->matrix, pDevice->rhs, incVsb, NIL(spREAL), NIL(spREAL));
    
    for ( eIndex = 1; eIndex <= pDevice->numElems; eIndex++ ) {
      pElem = pDevice->elements[ eIndex ];
      for ( index = 0; index <= 3; index++ ) {
	if ( pElem->evalNodes[ index ] ) {
	  pNode = pElem->pNodes[ index ];
	  if ( pNode->nodeType != CONTACT ) {
	    delPsi = incVsb[ pNode->psiEqn ] * delVsb;
	    solution[ pNode->psiEqn ] += delPsi;
	    if ( pElem->elemType == SEMICON
		&& (!OneCarrier || (OneCarrier == N_TYPE)) ) {
	      delN = incVsb[ pNode->nEqn ] * delVsb;
	      nConc = solution[ pNode->nEqn ];
	      newN = nConc + delN;
	      if ( newN <= 0.0 ) {
		solution[ pNode->nEqn ] = guessNewConc( nConc, delN );
	      }
	      else {
		solution[ pNode->nEqn ] = newN;
	      }
	    }
	    if ( pElem->elemType == SEMICON
		&& (!OneCarrier || (OneCarrier == P_TYPE)) ) {
	      delP = incVsb[ pNode->pEqn ] * delVsb;
	      pConc = solution[ pNode->pEqn ];
	      newP = pConc + delP;
	      if ( newP <= 0.0 ) {