Mercurial > forge
changeset 12117:d830ec7b1ca6 octave-forge
added circuit variables output to secs3d_coupled_circuit_newton.m
author | dvd7587 |
---|---|
date | Fri, 25 Oct 2013 10:33:29 +0000 |
parents | 88c10b8baba2 |
children | ecb65e433faf |
files | extra/secs1d/inst/demos/coupled_3d_diode.m extra/secs1d/inst/demos/coupled_3d_resistor_and_circuit.m extra/secs1d/inst/secs3d_coupled_circuit_newton.m |
diffstat | 3 files changed, 50 insertions(+), 34 deletions(-) [+] |
line wrap: on
line diff
--- a/extra/secs1d/inst/demos/coupled_3d_diode.m Wed Oct 23 17:35:22 2013 +0000 +++ b/extra/secs1d/inst/demos/coupled_3d_diode.m Fri Oct 25 10:33:29 2013 +0000 @@ -122,7 +122,7 @@ save -binary -z datafile_rlc_circuit.octbin.gz device material constants algorithm tspan Vin nin pin % A B C x r -[V, n, p, Fn, Fp, Jn, Jp, Itot, tout] = secs3d_coupled_circuit_newton ... +[V, n, p, F, Fn, Fp, Jn, Jp, Itot, tout] = secs3d_coupled_circuit_newton ... (device, material, constants, algorithm, Vin, nin, pin, tspan, @vbcs);
--- a/extra/secs1d/inst/demos/coupled_3d_resistor_and_circuit.m Wed Oct 23 17:35:22 2013 +0000 +++ b/extra/secs1d/inst/demos/coupled_3d_resistor_and_circuit.m Fri Oct 25 10:33:29 2013 +0000 @@ -1,3 +1,4 @@ +#!/usr/bin/octave %% clear all; close all; @@ -119,10 +120,11 @@ pause %% (pseudo)transient simulation - -save -binary -z datafile_rlc_circuit.octbin.gz device material constants algorithm tspan Vin nin pin % A B C x r - -[V, n, p, Fn, Fp, Jn, Jp, Itot, tout] = secs3d_coupled_circuit_newton ... +[~, ~, ~, ~, Fin, ~] = vbcs (0); +device.msh = bim3c_mesh_properties (device.msh); +save -binary -z datafile_rlc_circuit.octbin.gz device material constants algorithm tspan Vin nin pin Fin % A B C x r +return +[V, n, p, F, Fn, Fp, Jn, Jp, Itot, tout] = secs3d_coupled_circuit_newton ... (device, material, constants, algorithm, Vin, nin, pin, tspan, @vbcs);
--- a/extra/secs1d/inst/secs3d_coupled_circuit_newton.m Wed Oct 23 17:35:22 2013 +0000 +++ b/extra/secs1d/inst/secs3d_coupled_circuit_newton.m Fri Oct 25 10:33:29 2013 +0000 @@ -1,7 +1,13 @@ -function [V, n, p, Fn, Fp, Jn, Jp, Itot, tout] = ... +function [V, n, p, F, Fn, Fp, Jn, Jp, Itot, tout] = ... secs3d_coupled_circuit_newton ... (device, material, constants, algorithm, Vin, nin, pin, tspan, va) + + if (isfield (algorithm, 'plotsOff')); + plotsOff = algorithm.plotsOff; + else + plotsOff = 0; + endif tags = {'V', 'n', 'p', 'F'}; rejected = 0; @@ -21,9 +27,9 @@ this_dnodes = bim3c_unknowns_on_faces ( device.msh, device.contacts(iii)); if (! sum(this_dnodes)); - error((["circuit pin #%d is not connected", - " to the device. Check mesh bound" - "ary labeling."])'(:)', iii) + error((['circuit pin #%d is not connected', + ' to the device. Check mesh bound' + 'ary labeling.'])'(:)', iii) endif inodes = setdiff (inodes, this_dnodes); dnodes(iii, this_dnodes) = 1; @@ -64,19 +70,23 @@ dt, A, B, C, r, pins, indexing, dnodes, inodes); - jac = __secs3d_newton_jacobian__ (device, material, constants, - algorithm, V2, n2, p2, F2, - dt, A, B, C, r, pins, - indexing, dnodes, inodes); + jac = __secs3d_newton_jacobian__ (device, material, constants, + algorithm, V2, n2, p2, F2, + dt, A, B, C, r, pins, + indexing, dnodes, inodes); J = sparse(rows(jac{1, 1}), columns(jac{1, 1})); + for iii = 1:4 for jjj = 1:4 J += jac{iii, jjj}; end end + save -binary -z octave_jacobian_and_residual.octbin.gz J res + return + delta = - J \ res; dn = dp = dV = zeros(rows (n), 1); @@ -102,7 +112,7 @@ tk = min ([tkv, tkn, tkp]); if (tk <= 0) - error ("relaxation parameter too small, die!") + error ('relaxation parameter too small, die!') endif V2 += tk * dV; n2 += tk * dn; @@ -127,8 +137,8 @@ [incr0, whichone] = max ([incr0v, incr0n, incr0p, incr0F]); if (incr0 > algorithm.maxnpincr) - printf ("at time step %d, fixed point iteration %d, the ", tstep, in); - printf ("increment in %s has grown too large\n", tags{whichone}); + printf ('at time step %d, fixed point iteration %d, the ', tstep, in); + printf ('increment in %s has grown too large\n', tags{whichone}); reject = true; break; endif @@ -148,22 +158,24 @@ [incr1, whichone] = max ([incr1v, incr1n, incr1p, incr1F]); resnlin(in) = incr1; if (in > 3 && resnlin(in) > resnlin(in-3)) - printf ("at time step %d, fixed point iteration %d, ", tstep, in); - printf ("the Newton algorithm is diverging: "); - printf ("the increment in %s is not decreasing\n", tags{whichone}); + printf ('at time step %d, fixed point iteration %d, ', tstep, in); + printf ('the Newton algorithm is diverging: '); + printf ('the increment in %s is not decreasing\n', tags{whichone}); reject = true; break; endif - figure (1) - semilogy (1:in, resnlin(1:in),'bo-'); - xlim([1,15]); - ylim([5e-9,5e-2]); - drawnow + if (! plotsOff) + figure (1) + semilogy (1:in, resnlin(1:in),'bo-'); + xlim([1,15]); + ylim([5e-9,5e-2]); + drawnow + endif if (incr1 < algorithm.toll) - printf ("fixed point iteration %d, time step %d, ", in, tstep); - printf ("model time %g: convergence reached incr = %g ", t, incr1); + printf ('fixed point iteration %d, time step %d, ', in, tstep); + printf ('model time %g: convergence reached incr = %g ', t, incr1); break; endif @@ -175,10 +187,10 @@ t = tout (--tstep); dt /= 2; - printf ("reducing time step\n"); - printf ("\ttime step #%d, ", tstep); - printf ("model time %g s\n", t); - printf ("\tnew dt %g s\n", dt); + printf ('reducing time step\n'); + printf ('\ttime step #%d, ', tstep); + printf ('model time %g s\n', t); + printf ('\tnew dt %g s\n', dt); else @@ -212,12 +224,14 @@ endfor Itot(:, tstep) = (Jn_pins + Jp_pins + Jd_pins)(:); - figure (2) - plotyy (tout, Itot(2, :), tout, F(pins(2), :) - F(pins(1), :)); - drawnow + if (! plotsOff) + figure (2) + plotyy (tout, Itot(2, :), tout, F(pins(2), :) - F(pins(1), :)); + drawnow + endif dt *= .8 * sqrt (algorithm.maxnpincr / incr0); - printf ("\nestimate for next time step size: dt = %g \n", dt); + printf ('\nestimate for next time step size: dt = %g \n', dt); endif endwhile %% time step