Implementation of Source Concept in Matlab Miloslav Capek - - PowerPoint PPT Presentation

Implementation of Source Concept in Matlab

Miloslav ˇ Capek

Department of Electromagnetic Field CTU in Prague, Czech Republic miloslav.capek@fel.cvut.cz

Seminar at KTH Stockholm, Sweden January 19, 2017

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 1 / 31

Outline

1 Source Concept 2 Yet Another Challenge 3 Method of Moments 4 AToM 5 AToM Architecture 6 Useful Tools to Know 7 Conclusion In this talk: ◮ for the first talk in the series proceed to capek.elmag.org, ◮ elementary knowledge in Matlab is expected, ◮ be extremely careful when comparing different sources (papers): .

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 2 / 31

different notation!

Source Concept

Source Concept

Source Concept

Integral and variational methods Modal de- compositions Perspective topol-

• gy and

geometry HPC, algorithm efficiency Heuristic

• r convex
• ptimization

Sketch of main fields of the source concept.

Source concept. . . ◮ represents a radiator(s) completely by source currents J and M, ◮ infers all parameters from a source current, i.e. L (J , M).

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 3 / 31

Source Concept

What operators we know?

◮ assume operators in their matrix forms, i.e., L → L ⇒ J → I

• as far as meaningful, all algebraic operations are explicitly executable

Impedance op. Pr + 2jω (Wm − We) → 1 2IHZI

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 4 / 31

Source Concept

What operators we know?

◮ assume operators in their matrix forms, i.e., L → L ⇒ J → I

• as far as meaningful, all algebraic operations are explicitly executable

Impedance op. Pr + 2jω (Wm − We) → 1 2IHZI Stored energy op. (ka < 1) Wsto → 1 4IH ∂X ∂ω I = 1 4IHX′I

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 4 / 31

Source Concept

What operators we know?

◮ assume operators in their matrix forms, i.e., L → L ⇒ J → I

• as far as meaningful, all algebraic operations are explicitly executable

Impedance op. Pr + 2jω (Wm − We) → 1 2IHZI Stored energy op. (ka < 1) Wsto → 1 4IH ∂X ∂ω I = 1 4IHX′I Far-field op. Fι (r) → Fι (r) I Near-field op. E (r) → Ne (r) I, H (r) = Nm (r) I Directivity op. Dι (r) → 4πIHFH

ι (r) Fι (r) I

Z0IHRI = IHUι (r) I IHRI

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 4 / 31

Source Concept

What operators we know?

◮ assume operators in their matrix forms, i.e., L → L ⇒ J → I

• as far as meaningful, all algebraic operations are explicitly executable

Impedance op. Pr + 2jω (Wm − We) → 1 2IHZI Stored energy op. (ka < 1) Wsto → 1 4IH ∂X ∂ω I = 1 4IHX′I Far-field op. Fι (r) → Fι (r) I Near-field op. E (r) → Ne (r) I, H (r) = Nm (r) I Directivity op. Dι (r) → 4πIHFH

ι (r) Fι (r) I

Z0IHRI = IHUι (r) I IHRI Resistance op. PL → 1 2IHΣI Electric moment op. p → PI Magnetic moment op. m → MI

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 4 / 31

Source Concept

What are Operators Good For?

W J1 W J2 WB WA J(rA) J(rB) W0 Wmax Wfinal Wmax

VG1 VG2A VG2B

• 1. Determine lower/upper bounds,
• 2. estimate optimal trade-offs between

various parameters,

• 3. find proper feeding networks,
• 4. understand underlying physics,
• 5. try to make material synthesis,
• 6. try to modify shape of source region.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 5 / 31

Source Concept

Example: Characteristic Modes of PEC Fractal Plate

XIu = λuRIu

IFS fractal, 2nd iteration, PEC, ka = 0.5.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 6 / 31

Source Concept

Example: Characteristic Modes of PEC Fractal Plate

XIu = λuRIu

Current density J1 of the first characteristic mode. Radiation pattern for the first characteristic mode.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 6 / 31

Source Concept

Example: Spatial Decomposition of U-notched antenna

1 2 IHRI = 1 2

• IH

1

IH

2

R11 R12 R21 R22 I1 I2

• ka

0.6 0.8 1 1.2 1.4 Pr [W]

• 2
• 1

1 2 3

entire antenna arms only meanders only cross-terms Pr

(J1ÈJ2, J1ÈJ2)

Pr

(J1, J1)

Pr

(J2, J2)

2Pr

(J1, J2)

J1ÈJ2 J1 J2

Spatial (structural) decomposition of double U-notched antenna.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 7 / 31

Yet Another Challenge

Feeding Synthesis

◮ optimal currents are not compatible with realistic feeding

• check that V = ZIopt is dense vector full of non-zero entries

W

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 8 / 31

Yet Another Challenge

Feeding Synthesis

◮ optimal currents are not compatible with realistic feeding

• check that V = ZIopt is dense vector full of non-zero entries

◮ structure needs to be modified

• currents sub-optimal to Iopt
• heuristic optimization needed

W

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 8 / 31

Yet Another Challenge

Feeding Synthesis

◮ optimal currents are not compatible with realistic feeding

• check that V = ZIopt is dense vector full of non-zero entries

◮ structure needs to be modified

• currents sub-optimal to Iopt
• heuristic optimization needed

◮ shape modification resembles NP-hard problem

How much DOF we have?

W

N (unknowns) 28 52 120 ∞ possibilities unique solutions

Complexity of geometrical optimization for given voltage gap (red line) and N unknowns.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 8 / 31

Yet Another Challenge

Feeding Synthesis

◮ optimal currents are not compatible with realistic feeding

• check that V = ZIopt is dense vector full of non-zero entries

◮ structure needs to be modified

• currents sub-optimal to Iopt
• heuristic optimization needed

◮ shape modification resembles NP-hard problem

How much DOF we have?

W

N (unknowns) 28 52 120 ∞ possibilities 5.24 · 1029 1.39 · 1068 1.15 · 10199 ∞ unique solutions 2.68 · 108 4.50 · 1015 1.33 · 1036 ∞

Complexity of geometrical optimization for given voltage gap (red line) and N unknowns.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 8 / 31

Yet Another Challenge

Pixeling with Heuristic Optimization

Example: Shape optimization (pixeling) with ad hoc specified feeding placement, PEC plate L × L/2, ka = 0.3.

Antenna synthesis – how far can we go? ◮ On the present, only the heuristic optimization1. . .

• 1Y. Rahmat-Samii and E. Michielssen, Eds., Electromagnetic optimization by genetic algorithm. Wiley, 1999

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 9 / 31

Yet Another Challenge

Pixeling with Heuristic Optimization

Example: Shape optimization (pixeling) with ad hoc specified feeding placement, PEC plate L × L/2, ka = 0.3.

Antenna synthesis – how far can we go? ◮ On the present, only the heuristic optimization1. . . Computational time: 12116 s

Result of heuristic structural optimization using MOGA NSGAII from AToM-FOPS.

• 1Y. Rahmat-Samii and E. Michielssen, Eds., Electromagnetic optimization by genetic algorithm. Wiley, 1999

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 9 / 31

Yet Another Challenge

Pixeling with Heuristic Optimization

Example: Shape optimization (pixeling) with ad hoc specified feeding placement, PEC plate L × L/2, ka = 0.3.

Antenna synthesis – how far can we go? ◮ On the present, only the heuristic optimization1. . . Computational time: 12116 s

Result of heuristic structural optimization using MOGA NSGAII from AToM-FOPS.

Computational time: 1155 s

Result of deterministic in-house algorithm removing in each iteration the “worse” edge.

• 1Y. Rahmat-Samii and E. Michielssen, Eds., Electromagnetic optimization by genetic algorithm. Wiley, 1999

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 9 / 31

Yet Another Challenge

Pixeling with Heuristic Optimization

Example: Shape optimization (pixeling) with ad hoc specified feeding placement, PEC plate L × L/2, ka = 0.3.

Antenna synthesis – how far can we go? ◮ On the present, only the heuristic optimization1. . . Q (I) /QTM

Chu = 7.23

Resulting sub-optimal current approaching minimal value of quality factor Q.

Q (I) /QTM

Chu = 7.24

Resulting current given by in-house deterministic algorithm.

• 1Y. Rahmat-Samii and E. Michielssen, Eds., Electromagnetic optimization by genetic algorithm. Wiley, 1999

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 9 / 31

Yet Another Challenge

Structure of Solution Space

◮ all combinations for N = 28 edges (5.24 · 1029) calculated in Matlab

2

• 3 days on supercomputer, 2 resulting vectors + permutation table ≈ 55 GB

× × × × × × ×

2Work still in progress. For recent updates visit EuCAP 2017 presentation Excitation of Optimal and Suboptimal Currents. ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 10 / 31

Yet Another Challenge

Structure of Solution Space

◮ all combinations for N = 28 edges (5.24 · 1029) calculated in Matlab

2

• 3 days on supercomputer, 2 resulting vectors + permutation table ≈ 55 GB

280 300 320 340 360 380 400 7×104 6×104 5×104 4×104 3×104 2×104 1×104 quality factor Q number of solutions

best solution 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29

Structure of all suboptimal solution within 2 % tolerance to the best found candidate. Edge no. 18 is fed. 1 29 5 10 20 25 edge solution

fed edge removed edge retained edge best solution

1 25 75 Number of solutions in dependence on their quality factor Q. The best solution reaches Q (Ωopt) ≈ 292.

2Work still in progress. For recent updates visit EuCAP 2017 presentation Excitation of Optimal and Suboptimal Currents. ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 10 / 31

Yet Another Challenge

Example: Pareto Fronts For Excitation of Array #1

2 3 4 5 6 7 8 9 10 2 4 6 8 max {D} min {Q}

K = 4, B = 6 K = 6, B = 6 K = 8, B = 6 K = 8, B = 8 K = 8, B = 12

Vk = skAkejϕk

sk(1) ∈ {−1, 1} Ak(B) ∈ [0, 1] ϕk(B) ∈ [−π/2, π/2]

Multi-criteria optimization of K dipoles of length L = λ/2, placed side-by-side and fed with the voltage gaps into the middle. Separation distance is d = L/4, optimization genes for polarity sk, amplitude Ak and phase ϕk are composed of 1/B/B bits.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 11 / 31

Yet Another Challenge

Example: Pareto Fronts For Excitation of Array #2

Multi-criteria optimization of 4 dipoles of length L = λ/2, placed side-by-side and fed with the voltage gaps into the middle. Separation distance is d = L/4 (blue), d = L/8 (green), and d = L/16 (red), respectively. Optimization genes for polarity sk, amplitude Ak and phase ϕk of the driven voltage are composed of 1/8/8 bits.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 12 / 31

1 2 3 4 5 0.2 0.4 0.6 0.8 1 max {D} max {ηrad}

Yet Another Challenge

Key Functionality

◮ good mesh discretization is half of success ◮ all operators require extensive work with the mesh grid and its post-processing3 (MoM) ◮ most important operators (R, X, X′) are closely related to the method of moments4 ◮ various modal decompositions need tracking algorithm5 ◮ algorithms for single- or multi-criteria heuristic optimization6

• 4R. F. Harrington, Field computation by moment methods. Wiley – IEEE Press, 1993
• 3S. M. Rao, D. R. Wilton, and A. W. Glisson, “Electromagnetic scattering by surfaces of arbitrary shape”,

IEEE Trans. Antennas Propag., vol. 30, no. 3, pp. 409–418, 1982. doi: 10.1109/TAP.1982.1142818

• 5M. Capek, P. Hazdra, P. Hamouz, et al., “A method for tracking characteristic numbers and vectors”,

Prog.

• Electromagn. Res. B, vol. 33, pp. 115–134, 2011. doi: 10.2528/PIERB11060209
• 6K. Deb, Multi-objective optimization using evolutionary algorithms. Wiley, 2001

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 13 / 31

Method of Moments

MoM: Practical Aspects to Take Care Of

◮ discretization

• rectangles / triangles / quadrilateral elements
• in-house / FEKO export / Gmsh + MeshLab

◮ basis functions

• roof-tops / RWG / higher-order

◮ singularity treatment

• semi-numerical / numerical / . . .

◮ Gaussian quadrature

• higher-order / adaptive

◮ testing procedure

• Galerkin / non-Galerkin

◮ advanced features

• internal MoM normalization
• separation to “boxes”
• block vectorization
• . . .

11 5 7 8 8 4 4 6 2 3 12 3 10 2 7 3 5 6 9 z 1 1 y 2 1 x

Test case for MoM symmetries.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 14 / 31

Method of Moments

Useful Benchmark for Method of Moments

Spherical shell of radius a.

◮ synthetic test of quality of impedance matrices needed ◮ few fully separable bodies of finite extent ◮ CMs of spherical shell constitute excellent candidate

7

◮ characteristic eigenvalues λTE

n

= −yn (ka) jn (ka) , (1) λTM

n

= −(n + 1) yn (ka) − ka yn+1 (ka) (n + 1) jn (ka) − ka jn+1 (ka) , (2)

7See elmag.fel.cvut.cz/CMbenchmark. A paper proposing the benchmarks, their results and comparison of

state-of-the-art simulators will be published soon.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 15 / 31

x y z r ϕ ϑ

Method of Moments

Eigenvalues λn for CMs of Spherical Shell

0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 −20 −15 −10 −5 5 10 15 20 Wm < We TM modes ka λn

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 16 / 31

Method of Moments

Eigenvalues λn for CMs of Spherical Shell

0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 −20 −15 −10 −5 5 10 15 20 Wm < We TM modes TE modes Wm > We ka λn

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 16 / 31

Method of Moments

Example: Comparison of SW Packages

80 63 48 35 24 15 8 3 3 8 15 24 35 48 63 80 5 10 15 TM/TE mode order log10 |λn| TM modes TE modes exact FEKO (Zsym) CEM One (3/6) AToM (4) Makarov (Zsym)

Characteristic numbers for spherical shell of electrical size ka = 0.5, discretized into 500 triangles.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 17 / 31

Method of Moments

Example: Dependence on Order of Quadrature

24 35 48 63 80 99 6 8 10 12 14 TE mode order log10 |λn| Z = R + jX 24 35 48 63 80 99 TE mode order Zpos = ˆ Iξˆ IT + jX exact 1 2 4 7 10

Comparison of the quality of decomposition for original matrix and matrix post-processed so that no negative eigenvalues of RI = ξnI are present.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 18 / 31

Method of Moments

Example: Practical Indefiniteness of R

100 200 300 400 500 600 700 −15 −10 −5 mode index n log10 |ξn|

I II

RIn = ξnIn ξn > 0 ξn < 0 FEKO (original data) FEKO (symmetrized) AToM (1st order) AToM (8th order)

Eigenvalues for decomposition of real part of impedance matrix for spherical shell, discretized into 500 triangles.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 19 / 31

AToM

AToM: Antenna Toolbox For Matlab

“Antenna source concept” – New approach to antenna design.

AToM (Antenna Toolbox for Matlab) third-party product developed by CTU in Prague and BUT.

Logo of the AToM project.

The main idea behind the AToM toolbox is to develop new package that will be able to: ◮ utilize the source concept features ◮ handle with data from third party software ◮ accept other codes from the community ◮ make it possible the fast-prototyping of advanced antenna designs

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 20 / 31

AToM Architecture

Matlab-like Conception

User AToM: GUI High-level functions Low-level functions Skilled user Expert user MATLAB

Scheme of AToM – (almost) completely written in OOP.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 21 / 31

AToM Architecture

Why Matlab?

Pros ◮ high-definition language

• excellent for fast-prototyping
• many built-in functions are

embedded

◮ new functionality can easily be published8 ◮ maybe other. . . Cons ◮ still not as fast as e.g. C

• and to be efficient, Matlab

needs very good programming skill

◮ not open-source ◮ to make standalone application is a nightmare ◮ maybe other. . . ◮ What is your opinion?

8www.mathworks.com/matlabcentral/fileexchange ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 22 / 31

Useful Tools to Know

Features of Matlab (>R2015b)

◮ Run-time Type Analysis (>Matlab 6.5)

• data types in m-file are noticed during the first run

◮ Just-In-Time-Accelerator9

• parts of code that satisfy certain conditions are precompiled
• runs always (>Matlab 2016a)

◮ profiling via profile

• JIT however deactivated during the profile measurement
• almost impossible to do good job without it

◮ unit-test framework (>2014b) ◮ Source Control Integration

• GIT, SVN
• Jenkins can be utilized

9See http://www.ee.columbia.edu/ marios/matlab/accel matlab.pdf. ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 23 / 31

Useful Tools to Know

Agile/Team Development

02/2015 01/2016 01/2017 20000 40000 60000 80000 100000 120000 time number of lines lines of code 02/2015 01/2016 01/2017 1000 2000 3000 time number of m-files and m-functions m-files m-functions

Some statistics of AToM project over time (Matlab+GIT+Jenkins).

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 24 / 31

Useful Tools to Know

Agile/Team Development

02/2015 01/2016 01/2017 20000 40000 60000 80000 100000 120000 time number of lines lines of code 02/2015 01/2016 01/2017 1000 2000 3000 time number of m-files and m-functions m-files m-functions

Some statistics of AToM project over time (Matlab+GIT+Jenkins).

directories 136 packages 83 classes 121 m-files 1659 functions 2764 unitTests 1176 lines of code 123280 comments 12441

Valid data on 19/01/2017, 00:23 AM.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 24 / 31

Useful Tools to Know

Object-Oriented Programming (OOP)

◮ Object-oriented programming (>R2008)

• surprisingly rich OOP (all classical OOP patterns feasible)
• starts to be integrated everywhere (see, e.g., new graphics in >R2014a)
• numeric × handle classes
• multiple inheritance (no interfaces etc.)
• yED for UML

◮ sustainable code development

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 25 / 31

Useful Tools to Know

OOP & vectorization

1 % ... 2 % check whos weight is even number 3 function val = isWeightEven(objs) 4 val = mod([objs.weight],2) == 0; 5 end 6 7 % return "object" who is the oldest one... 8 function oldest obj = whoIsTheOldest(objs) 9 allAges = [objs.age]; % for acceleration purposes 10

• ldest obj = objs(allAges == max(allAges));

11 end 12 13 % increase age of all objects 14 function objs = increaseAge(objs, incr age) 15 [objs.age] = indexing.listEntries([objs.age] + incr age); 16 end 17 % ... 1 % increase age (modification) − FOR approach 2 for thisObj = 1:N(thisN) 3 ppl(thisObj).increaseAge(10); 4 end 1 % increase age (modification) − vectorized approach 2 ppl.increaseAge(10);

OOP and vectorization (highest level of abstraction in Matlab).

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 26 / 31

Useful Tools to Know

Speed of Code

101 102 103 104 10−4 10−3 10−2 10−1 100 Size of date set Evaluation time [s] (warm-up done) R2011a 101 102 103 104 Size of date set R2014b 101 102 103 104 Size of date set R2016b for creation comparison indexing modification vect. creation comparison indexing modification

Matlab implementation of OOP: comparison of for cycles with vectorized OOP.

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 27 / 31

Useful Tools to Know

New Features

◮ from R2017a the functions can be inside scripts ◮ OOP improvements (e.g., abstract classes can have constructors) ◮ full OOP graphics

• figure contains all the data
• number and position of markers can be specified
• many small but useful changes

◮ big data (tallArrays, saveValue) ◮ implicitly applied repmat function ◮ Matlab Live Scripts ◮ uifigure and ui-objects

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 28 / 31

Conclusion

About CTU in Prague

Established in 1707 as the first non-military university in Central Europe. ◮ from 12 students in 1707 to more than 20 000 students in 2016

CTU, Faculty of Electrical Engineering (one of eight faculties).

You are welcome to visit us in Prague!

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 29 / 31

Conclusion

Conclusion

◮ nowadays, good physicists/engineers should be good programmers as well

• from implementation of theory to controlling the measuring devices

effectively

• Matlab/Python+Mathematica, LaTeX (TikZ)

◮ knowledge in EM simulators is essential to check our calculations and implementations (simulations are like measurements in the old days) ◮ understanding of numerical methods (at least one) is advantage

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 30 / 31

Questions?

For complete PDF presentation see

capek.elmag.org

Miloslav ˇ Capek miloslav.capek@fel.cvut.cz

• 19. 1. 2017, v1.3

ˇ Capek, M., CTU in Prague Implementation of Source Concept in Matlab 31 / 31