From edaf2463e8ae05b00b659ecbe46b21ab0e99c712 Mon Sep 17 00:00:00 2001 From: Martin Holters Date: Wed, 4 Aug 2021 08:49:05 +0200 Subject: [PATCH] Fix Jiles-Atherton-based inductor/transformer model And add tests to ensure they don't bit-rot again. --- src/elements.jl | 17 ++++++++-------- test/runtests.jl | 51 ++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 60 insertions(+), 8 deletions(-) diff --git a/src/elements.jl b/src/elements.jl index 7781d52f..54100d71 100644 --- a/src/elements.jl +++ b/src/elements.jl @@ -1,4 +1,4 @@ -# Copyright 2015, 2016, 2017, 2018, 2019, 2020 Martin Holters +# Copyright 2015, 2016, 2017, 2018, 2019, 2020, 2021 Martin Holters # See accompanying license file. export resistor, potentiometer, capacitor, inductor, transformer, @@ -57,7 +57,8 @@ Henri. The coupling can either be specified using `coupling_coefficient` (0 is not coupled, 1 is closely coupled) or by `mutual_coupling`, the mutual inductance in Henri, where the latter takes precedence if both are given. -Pins: `1` and `2` for primary winding, `3` and `4` for secondary winding +Pins: `primary1` and `primary2` for primary winding, `secondary1` and +`secondary2` for secondary winding """ transformer(l1, l2; coupling_coefficient=1, mutual_coupling=coupling_coefficient*sqrt(l1*l2)) = @@ -91,7 +92,7 @@ of nonlinear inductor circuits,” IEEE Trans. Magn., vol. 30, no. 5, pp. 2795–2801, 1994, where the definition of `c` is taken from the latter. The ACME implementation is discussed in [M. Holters, U. Zölzer, "Circuit Simulation with Inductors and Transformers Based on the Jiles-Atherton Model of -Magnetization"](http://ant-s4.unibw-hamburg.de/paper-archive/2016/dafxpapers/08-DAFx-16_paper_10-PN.pdf). +Magnetization"](http://dafx.de/paper-archive/2016/dafxpapers/08-DAFx-16_paper_10-PN.pdf). Pins: `1` and `2` for primary winding, `3` and `4` for secondary winding, and so on @@ -120,13 +121,13 @@ function transformer(::Type{Val{:JA}}; D=2.4e-2, A=4.54e-5, ns=[], J_1_2 = (1e-4/Ms) * -(1-c)^2*k * δM*δ/den^2 * q[3] J_1_3 = (1e-4/Ms) * ((1-c) * δM*(Man-q[2])/den + (c*Ms/a)*Ld_q1) J_1_4 = (1e-4/Ms) * ((c * Ms/a * α)*Ld_q1 - 1) - return (res, @SMatrix [J_1_1 J_1_2; J_2_1 J_2_2]) + return (res, @SMatrix [J_1_1 J_1_2 J_1_3 J_1_4]) end - Element(mv=[speye(length(ns)); spzeros(5, length(ns))], + Element(mv=Matrix{Int}(I, length(ns)+5, length(ns)), mi=[spzeros(length(ns), length(ns)); ns'; spzeros(4, length(ns))], mx=[spzeros(length(ns), 2); -π*D 0; -1/a -α/a; 0 -1; 0 0; 0 0], mxd=[-μ0*A*ns -μ0*ns*A; 0 0; 0 0; 0 0; -1 0; 0 -1], - mq=[zeros(length(ns)+1,4); eye(4)], nonlinear_eq = nonlinear_eq) + mq=[zeros(length(ns)+1,4); Matrix(I, 4, 4)], nonlinear_eq = nonlinear_eq) end @@ -148,14 +149,14 @@ parameters are set using named arguments: | `k` | amount of hysteresis (in Ampere-per-meter) | | `Ms` | saturation magnetization (in Ampere-per-meter) | -A detailed discussion of the paramters can be found in D. C. Jiles and D. L. +A detailed discussion of the parameters can be found in D. C. Jiles and D. L. Atherton, “Theory of ferromagnetic hysteresis,” J. Magn. Magn. Mater., vol. 61, no. 1–2, pp. 48–60, Sep. 1986 and J. H. B. Deane, “Modeling the dynamics of nonlinear inductor circuits,” IEEE Trans. Magn., vol. 30, no. 5, pp. 2795–2801, 1994, where the definition of `c` is taken from the latter. The ACME implementation is discussed in [M. Holters, U. Zölzer, "Circuit Simulation with Inductors and Transformers Based on the Jiles-Atherton Model of -Magnetization"](http://ant-s4.unibw-hamburg.de/paper-archive/2016/dafxpapers/08-DAFx-16_paper_10-PN.pdf). +Magnetization"](http://dafx.de/paper-archive/2016/dafxpapers/08-DAFx-16_paper_10-PN.pdf). Pins: `1`, `2` """ diff --git a/test/runtests.jl b/test/runtests.jl index eda5a08f..c0182a6a 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -412,6 +412,57 @@ end @test run!(model, zeros(0,1)) ≈ [10/100000] end +@testset "Jiles-Atherton inductor/transformer" begin + @testset "Jiles-Atherton inductor" begin + circ = @circuit begin + Jin = voltagesource() + Jout1 = currentprobe(), [+] ↔ Jin[+] + Jout2 = currentprobe(), [+] ↔ Jin[+] + L_JA = inductor(Val{:JA}), [1] ↔ Jout1[-], [2] ↔ Jin[-] + # default parameters approximately linearize as 174mH + L_lin = inductor(174e-3), [1] ↔ Jout2[-], [2] ↔ Jin[-] + end + model = DiscreteModel(circ, 1//44100) + # starting non-magnetized, the inductor first behaves sub-linear + y = run!(model, fill(0.1, 1, 750)) + @test isapprox(y[1,1:9], y[2,1:9], rtol=1e-2) # almost linear at first + @test all(y[1,:] .< y[2,:]) + # until it reaches saturation, where it becomes super-linear + run!(model, fill(0.1, 1, 500)) + y = run!(model, fill(0.1, 1, 750)) + @test all(y[1,:] .> y[2,:]) + # due to hysteresis, applying the negative voltage for the same time + # drives current below zero + y = run!(model, fill(-0.1, 1, 2000)) + @test y[1,end] < -2e-3 + # with voltage at zero (i.e. shorted), the current stays constant + y = run!(model, zeros(1, 1000)) + @test y[1,1] < -2e-3 && all(y[:,1] .≈ y) + end + @testset "Jiles-Atherton transformer" begin + circ = @circuit begin + Jin = voltagesource() + R1 = resistor(10), [1] ↔ Jin[+] + R2 = resistor(10), [1] ↔ Jin[+] + T_JA = transformer(Val{:JA}; ns = [10, 100]), [1] ↔ R1[2], [2] ↔ Jin[-] + # approximate small-scale equivalent + T_lin = transformer(330e-6, 33e-3), [primary1] ↔ R2[2], [primary2] ↔ Jin[-] + Jout1 = voltageprobe(;gp=1e-3), [+] ↔ T_JA[3], [-] ↔ T_JA[4] + Jout2 = voltageprobe(;gp=1e-3), [+] ↔ T_lin[secondary1], [-] ↔ T_lin[secondary2] + end + model = DiscreteModel(circ, 1//44100) + u = [sin(2pi*1000/44100*n) for n in (0:499)]' + # almost linear for small input + y = run!(model, 0.001*u)[:,200:end] + @test isapprox(y[1,:], y[2,:]; rtol = 1.2e-3) + y = run!(model, 0.002*u)[:,200:end] + @test isapprox(y[1,:], y[2,:]; rtol = 1.2e-3) + # not at all linear for large input + y = run!(model, 10*u)[200:end] + @test !isapprox(y[1,:], y[2,:]; rtol = 0.75) + end +end + @testset "BJT" begin isc=1e-6 ise=2e-6