1 function [Fx,Fy]=Ufvsgcurrent3(mesh,u,alpha,gamma,eta,beta);
6 ## @deftypefn {Function File} {[@var{jx},@var{jy}]} = Ufvsgcurrent3 @
7 ## (@var{mesh}, @var{u}, @var{alpha}, @var{gamma}, @var{eta}, @var{beta});
10 ## Builds the Scharfetter-Gummel approximation of the vector
15 ## $ \vect{J}(u) = \alpha \gamma (\eta\vect{\nabla}u-\vect{beta}u) $
19 ## J(@var{u}) = @var{alpha}* @var{gamma} * (@var{eta} * grad @var{u} - @var{beta} * @var{u}))
24 ## @item @var{alpha} is an element-wise constant scalar function
25 ## @item @var{eta}, @var{u}, @var{gamma} are piecewise linear
26 ## conforming scalar functions
27 ## @item @var{beta} is an element-wise constant vector function
30 ## J(@var{u}) is an element-wise constant vector function
32 ## Instead of passing the vector field @var{beta} directly
33 ## one can pass a piecewise linear conforming scalar function
34 ## @var{phi} as the last input. In such case @var{beta} = grad @var{phi}
35 ## is assumed. If @var{phi} is a single scalar value @var{beta}
36 ## is assumed to be 0 in the whole domain.
38 ## @seealso{Uscharfettergummel3}
41 %% This file is part of
43 %% SECS2D - A 2-D Drift--Diffusion Semiconductor Device Simulator
44 %% -------------------------------------------------------------------
45 %% Copyright (C) 2004-2006 Carlo de Falco
49 %% SECS2D is free software; you can redistribute it and/or modify
50 %% it under the terms of the GNU General Public License as published by
51 %% the Free Software Foundation; either version 2 of the License, or
52 %% (at your option) any later version.
54 %% SECS2D is distributed in the hope that it will be useful,
55 %% but WITHOUT ANY WARRANTY; without even the implied warranty of
56 %% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
57 %% GNU General Public License for more details.
59 %% Youx should have received a copy of the GNU General Public License
60 %% along with SECS2D; If not, see <http://www.gnu.org/licenses/>.
62 Nelem = columns(mesh.t);
63 Nnodes = columns(mesh.p);
65 uloc = u(mesh.t(1:3,:));
67 shgx = reshape(mesh.shg(1,:,:),3,Nelem);
68 shgy = reshape(mesh.shg(2,:,:),3,Nelem);
70 x = reshape(mesh.p(1,mesh.t(1:3,:)),3,[]);
71 dx = [ (x(3,:)-x(2,:)) ;
75 y = reshape(mesh.p(2,mesh.t(1:3,:)),3,[]);
76 dy = [ (y(3,:)-y(2,:)) ;
80 ##ei = sqrt( dx.^2 + dy.^2 );
84 elseif all(size(beta)==[2,Nelem])
85 v23 = beta(1,:) .* dx(1,:) + beta(2,:) .* dy(1,:);
86 v31 = beta(1,:) .* dx(2,:) + beta(2,:) .* dy(2,:);
87 v12 = beta(1,:) .* dx(3,:) + beta(2,:) .* dy(3,:);
88 elseif all(size(beta)==[Nnodes,1])
89 betaloc = beta(mesh.t(1:3,:));
90 v23 = betaloc(3,:)-betaloc(2,:);
91 v31 = betaloc(1,:)-betaloc(3,:);
92 v12 = betaloc(2,:)-betaloc(1,:);
95 etaloc = eta(mesh.t(1:3,:));
97 eta23 = etaloc(3,:)-etaloc(2,:);
98 eta31 = etaloc(1,:)-etaloc(3,:);
99 eta12 = etaloc(2,:)-etaloc(1,:);
101 etalocm1 = Utemplogm(etaloc(2,:),etaloc(3,:));
102 etalocm2 = Utemplogm(etaloc(3,:),etaloc(1,:));
103 etalocm3 = Utemplogm(etaloc(1,:),etaloc(2,:));
105 gammaloc = gamma(mesh.t(1:3,:));
106 geloc = gammaloc.*etaloc;
108 gelocm1 = Utemplogm(geloc(2,:),geloc(3,:));
109 gelocm2 = Utemplogm(geloc(3,:),geloc(1,:));
110 gelocm3 = Utemplogm(geloc(1,:),geloc(2,:));
112 [bp23,bm23] = Ubern( (v23 - eta23)./etalocm1);
113 [bp31,bm31] = Ubern( (v31 - eta31)./etalocm2);
114 [bp12,bm12] = Ubern( (v12 - eta12)./etalocm3);
116 gfigfj = [ shgx(3,:) .* shgx(2,:) + shgy(3,:) .* shgy(2,:) ;
117 shgx(1,:) .* shgx(3,:) + shgy(1,:) .* shgy(3,:) ;
118 shgx(2,:) .* shgx(1,:) + shgy(2,:) .* shgy(1,:) ];
120 Fx = - alpha' .* ( gelocm1 .* etalocm1 .* dx(1,:) .* ...
122 ( bp23 .* uloc(3,:)./etaloc(3,:) -...
123 bm23 .* uloc(2,:)./etaloc(2,:)) +... %% 1
124 gelocm2 .* etalocm2 .* dx(2,:) .* ...
126 (bp31 .* uloc(1,:)./etaloc(1,:) -...
127 bm31 .* uloc(3,:)./etaloc(3,:)) +... %% 2
128 gelocm3 .* etalocm3 .* dx(3,:) .* ...
130 (bp12 .* uloc(2,:)./etaloc(2,:) -...
131 bm12 .* uloc(1,:)./etaloc(1,:)) ... %% 3
134 Fy = - alpha' .* ( gelocm1 .* etalocm1 .* dy(1,:) .* ...
136 ( bp23 .* uloc(3,:)./etaloc(3,:) -...
137 bm23 .* uloc(2,:)./etaloc(2,:)) +... %% 1
138 gelocm2 .* etalocm2 .* dy(2,:) .* ...
140 (bp31 .* uloc(1,:)./etaloc(1,:) -...
141 bm31 .* uloc(3,:)./etaloc(3,:)) +... %% 2
142 gelocm3 .* etalocm3 .* dy(3,:) .* ...
144 (bp12 .* uloc(2,:)./etaloc(2,:) -...
145 bm12 .* uloc(1,:)./etaloc(1,:)) ... %% 3