octave_packages/ocs-0.1.3/sbn/Mpdesympnjunct.m

   1 ## Copyright (C) 2006,2007,2008  Carlo de Falco
   2 ##
   3 ## This file is part of:
   4 ## OCS - A Circuit Simulator for Octave
   5 ##
   6 ## OCS is free software; you can redistribute it and/or modify
   7 ## it under the terms of the GNU General Public License as published by
   8 ## the Free Software Foundation.
   9 ##
  10 ## This program is distributed in the hope that it will be useful,
  11 ## but WITHOUT ANY WARRANTY; without even the implied warranty of
  12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  13 ## GNU General Public License for more details.
  14 ##
  15 ## You should have received a copy of the GNU General Public License
  16 ## along with this program (see the file LICENSE); if not,
  17 ## see <http://www.gnu.org/licenses/>.
  18 ##
  19 ## author: Carlo de Falco <cdf _AT_ users.sourceforge.net>
  20
  21 ## -*- texinfo -*-
  22 ##
  23 ## @deftypefn{Function File} @
  24 ## {[@var{j},@var{g}]=}Mpdesympnjunct@
  25 ## (@var{len}, @var{Dope}, @var{va}, @
  26 ## @var{Area}, @var{Nnodes}, @var{toll}, @var{maxit}, @var{ptoll}, @var{pmaxit})
  27 ##
  28 ## INTERNAL FUNCTION:
  29 ##
  30 ## NOT SUPPOSED TO BE CALLED DIRECTLY BY USERS
  31 ## @end deftypefn
  32
  33 function [j,g] = Mpdesympnjunct (len,Dope,va,Area,Nnodes,toll,maxit,ptoll,pmaxit)
  34
  35   constants
  36
  37   x = linspace(0,len,Nnodes)';
  38   xm = mean(x);
  39   Nd=Dope;
  40   Na=Dope;
  41   nn = (Nd + sqrt(Nd^2+4*ni^2))/2;
  42   pp = (Na + sqrt(Na^2+4*ni^2))/2;
  43   xn = xm+len/10;
  44   xp = xm-len/10;
  45   D = Nd * (x>xm) - Na * (x<=xm);
  46
  47   ## Scaling coefficients
  48   xs        = len;
  49   ns        = norm(D,inf);
  50   Vs        = Vth;
  51   us        = un;
  52   Js        = (Area*us*Vs*q*ns/xs);
  53   xin       = x/xs;
  54   gs        = Js/Vs;
  55
  56   n = nn * (x>=xn) + (ni^2)/pp * (x<xn);
  57   p = (ni^2)/nn * (x>xp) + pp * (x<=xp);
  58
  59   Fn = va*(x<=xm);
  60   Fp = Fn;
  61
  62   V  = (Fn - Vth * log(p/ni));
  63
  64   ## Scaling
  65   idata.D    = D/ns;
  66   idata.un   = un/us;
  67   idata.up   = up/us;
  68   idata.tn   = inf;
  69   idata.tp   = inf;
  70   idata.l2   = (Vs*esi)/(q*ns*xs^2);
  71   idata.nis  = ni/ns;
  72   idata.n   = n/ns;
  73   idata.p   = p/ns;
  74   idata.V   = V/Vs;
  75   idata.Fn  = (Fn - Vs * log(ni/ns))/Vs;
  76   idata.Fp  = (Fp + Vs * log(ni/ns))/Vs;
  77
  78
  79   ## Solution of DD system
  80   ## Algorithm parameters
  81
  82   sinodes = [1:length(x)];
  83
  84   [idata,it,res] = DDGgummelmap (xin,idata,toll,maxit/2,ptoll,pmaxit,0);
  85   [odata,it,res] = DDNnewtonmap (xin,idata,toll,maxit/2,0);
  86
  87   DV = diff(odata.V);
  88   h  = xin(2)-xin(1);
  89
  90   Bp = Ubernoulli(DV,1);
  91   Bm = Ubernoulli(DV,0);
  92
  93   Jn = -odata.un * (odata.n(2:end).*Bp-odata.n(1:end-1).*Bm)/h;
  94   Jp =  odata.up * (odata.p(2:end).*Bm-odata.p(1:end-1).*Bp)/h;
  95
  96   coeff = idata.un.*Umediaarmonica(odata.n);
  97   L =- Ucomplap (xin,Nnodes,[],[],coeff);
  98   Jn1 = L*odata.Fn;
  99   fprintf(1,"jn1=%g\n",Jn1(1))
 100   fprintf(1,"jn=%g\n",Jn(1))
 101
 102   C11 = L(1,1);
 103   C1I = L(1,2:end-1);
 104   CII = L(2:end-1,2:end-1);
 105   Gn   = C11 - C1I*(CII\C1I');
 106   Gn   = Gn - coeff(1)*(odata.Fn(2)-odata.Fn(1))/h
 107
 108   coeff = idata.up.*Umediaarmonica(odata.p);
 109   L =- Ucomplap (xin,Nnodes,[],[],coeff);
 110   Jp1 = L*odata.Fp;
 111   fprintf(1,"jp1=%g\n",Jp1(1))
 112   fprintf(1,"jp=%g\n",Jp(1))
 113
 114   C11 = L(1,1);
 115   C1I = L(1,2:end-1);
 116   CII = L(2:end-1,2:end-1);
 117   Gp   = C11 - C1I*(CII\C1I');
 118   Gp   = Gp - coeff(1)*(odata.Fp(2)-odata.Fp(1))/h
 119
 120
 121   ## Descaling
 122   j= -(Jp(1)+Jn(1))*Js
 123   g= gs*(Gn+Gp)
 124
 125 endfunction