Programming Fundamentals for Engineers

Chapter 13
Numerical and Data Programming for Electrical Engineering

Bachelors Program in Electrical Engineering
Programming Fundamentals for Engineers   •   Bachelor

Numerical and Data Programming for Electrical Engineering

Introduction

Electrical engineers routinely convert physical measurements into numerical models: voltages and currents become arrays, circuit equations become matrices, sensor samples become time series, and experimental uncertainty becomes statistical evidence. Programming makes these tasks repeatable, but numerical software is not automatically correct. A program may execute without error while producing an inaccurate or physically impossible result because of poor conditioning, inappropriate sampling, unit inconsistency, or floating-point round-off.

This chapter develops practical numerical programming skills for electrical engineering. Examples use Python-like notation with NumPy-style operations, but the mathematical ideas apply to MATLAB, Julia, C++, and other scientific-computing environments. The emphasis is on selecting suitable algorithms, checking assumptions, interpreting results, and communicating engineering conclusions.

Learning Objectives

After studying this chapter, you should be able to:

  1. Represent vectors, matrices, and complex quantities in engineering programs.
  2. Perform matrix operations and solve linear systems reliably.
  3. Calculate with complex impedance, phasors, and frequency response.
  4. Approximate integrals and find roots numerically.
  5. Process sampled signals and measurement data.
  6. Produce informative plots and statistical summaries.
  7. Evaluate floating-point error, conditioning, and numerical stability.
  8. Design reproducible and professionally defensible numerical workflows.

1. Numerical Representation in Engineering Programs

A numerical algorithm operates on finite representations of quantities that may be continuous or exact in theory. A measurement such as 3.30 V is stored as a binary floating-point approximation, not as an infinitely precise real number. Similarly, a continuous signal is represented by a finite sequence of samples.

1.1 Scalars, vectors, and matrices

A scalar represents one quantity. A vector represents an ordered collection, such as voltage samples:

v = [v0, v1, …, vN−1]T

Here, N is the number of samples and the superscript T denotes transpose. A matrix is a rectangular array. In circuit analysis, the nodal equation is commonly written as:

A x = b

where A is the coefficient matrix, x is the vector of unknown node voltages or currents, and b is the source vector. The units of rows and columns must be interpreted consistently. For example, a conductance matrix has units of siemens, a voltage vector has volts, and their product has amperes.

In Python-like code:

import numpy as np
A = np.array([[10.0, -2.0], [-2.0, 5.0]])
b = np.array([8.0, 3.0])
x = np.linalg.solve(A, b)

The preferred operation is solve(A, b), not inv(A) @ b. Explicit matrix inversion usually requires more operations and can amplify numerical error.

1.2 Elementwise and matrix operations

A frequent programming error is confusing elementwise multiplication with matrix multiplication. If u and v are arrays of equal shape, u * v multiplies corresponding elements. The matrix product is written A @ B, and its dimensions must satisfy the inner-dimension rule. If A is m × n and B is n × p, then AB is m × p.

The dot product of two real vectors is:

u · v = ∑k=0N−1 uk vk

For complex vectors, engineering software may require a conjugate inner product:

uᴴ v = ∑k=0N−1 uk* vk

where H denotes conjugate transpose and * denotes complex conjugation.

1.3 Matrix properties relevant to engineering

A square matrix A is nonsingular if its determinant is nonzero, equivalently if its inverse exists. However, determinant calculation is not a reliable general-purpose test of numerical solvability. A matrix can have a nonzero determinant but still be badly conditioned.

The condition number in a chosen norm is:

κ(A) = ||A|| ||A⁻¹||

For a nonsingular matrix, κ(A) ≥ 1. A large condition number indicates that small changes in input data may produce large changes in the solution. The condition number depends on the norm and on scaling. It is a sensitivity measure, not necessarily a direct measure of the error in one particular computation.

2. Complex Numbers and Phasor Calculations

Sinusoidal steady-state analysis represents a sinusoidal quantity by a complex phasor. Using RMS phasors, a voltage can be written as:

V = |V| ejθ = |V|(cos θ + j sin θ)

where |V| is RMS magnitude, θ is phase angle, and j is the engineering symbol for √−1. The corresponding time-domain voltage is:

v(t) = √2 |V| cos(ωt + θ)

This relation assumes a sinusoidal steady-state signal and RMS phasor magnitude. If peak phasors are used instead, the factor √2 is omitted.

For impedance:

ZR = R
ZL = jωL
ZC = 1/(jωC)

Here R is resistance in ohms, L is inductance in henries, C is capacitance in farads, and ω = 2πf is angular frequency in radians per second. These expressions assume linear components and sinusoidal steady state at frequency f.

The phasor form of Ohm’s law is:

V = ZI

For a complex number z = a + jb:

|z| = √(a² + b²)
arg(z) = atan2(b, a)

The two-argument arctangent is preferred because it identifies the correct quadrant.

Worked Example 1: Series RLC circuit

Consider a 10 Ω resistor, a 20 mH inductor, and a 100 μF capacitor in series at f = 50 Hz, supplied by 20 V RMS. First calculate:

ω = 2πf = 314.159 rad/s

The reactances are:

XL = ωL = 314.159 × 0.020 = 6.283 Ω
XC = 1/(ωC) = 1/(314.159 × 0.0001) = 31.831 Ω

Therefore:

Z = R + j(XL − XC) = 10 − j25.548 Ω

Its magnitude is approximately 27.435 Ω, so the RMS current magnitude is:

|I| = |V|/|Z| = 20/27.435 ≈ 0.729 A

The negative imaginary part indicates a capacitive net reactance, and the current leads the source voltage by approximately 68.6°. A program should preserve complex data types:

f = 50.0
w = 2*np.pi*f
Z = 10 + 1j*(w*0.020 - 1/(w*100e-6))
I = 20/Z

A common mistake is to use degrees in functions that expect radians. Convert explicitly when necessary.

3. Solving Linear Systems and Matrix Models

Many engineering models produce simultaneous linear equations. Examples include nodal analysis, mesh analysis, state-space models, least-squares estimation, and finite-difference approximations.

For a 2 × 2 system:

a11x1 + a12x2 = b1
a21x1 + a22x2 = b2

A numerical library normally applies factorization methods such as LU decomposition rather than symbolic elimination. After calculating x, verify the residual:

r = b − A x

A small residual is necessary but not always sufficient for an accurate solution when A is ill-conditioned. A useful relative residual is:

ρ = ||b − Ax|| / (||A|| ||x|| + ||b||)

The denominator prevents misleadingly large or small interpretations caused by scale. The norm must be specified; the Euclidean norm is common.

Least-squares estimation

When there are more measurements than unknown parameters, the system Ax ≈ b is usually solved by minimizing:

minx ||Ax − b||₂²

If A has full column rank, the normal-equation solution is:

x = (AᵀA)⁻¹ Aᵀb

This formula is mathematically valid under the stated rank condition, but directly forming AᵀA can worsen conditioning. QR factorization or singular-value decomposition is generally preferred in software. For complex data, transpose is replaced by conjugate transpose Aᴴ.

4. Numerical Integration and Root Finding

4.1 Numerical integration

An integral represents accumulated quantity. For example, energy delivered by a time-varying power signal is:

E = ∫t0t1 p(t) dt

If p is in watts and t is in seconds, E is in joules. With sampled values p_k at equally spaced times separated by Δt, the trapezoidal rule is:

∫t0tN p(t)dt ≈ Δt[½p0 + ∑k=1N−1pk + ½pN]

For a sufficiently smooth function, the composite trapezoidal error generally decreases proportionally to Δt² as the sampling interval decreases. The method is not automatically accurate for discontinuities or rapidly varying signals.

In Python-like code:

energy = np.trapezoid(power, time)

The time array should be monotonic, and nonuniform sampling must be handled by supplying the actual time values.

4.2 Root finding

A root of f(x) satisfies:

f(x) = 0

Roots may represent a circuit operating point, a resonant frequency, or the time at which a signal crosses a threshold.

The bisection method requires a continuous function on [a, b] with:

f(a)f(b) < 0

This sign change guarantees at least one root in the interval by the intermediate value theorem. At each step, c = (a + b)/2 is evaluated, and the half-interval retaining a sign change is selected. After n iterations, the interval width is (b − a)/2ⁿ.

Newton’s method uses:

xk+1 = xk − f(xk)/f′(xk)

It can converge rapidly near a simple root, but requires a derivative or derivative approximation and may diverge from a poor initial guess or near a point where f′(x) is small. Bracketing methods are usually safer for production engineering software.

5. Signal and Measurement Processing

5.1 Sampling and units

A sampled signal is represented as x[n] = x(nT_s), where T_s is the sampling interval and f_s = 1/T_s is the sampling frequency. To avoid aliasing for a strictly band-limited signal with highest frequency f_max, the sampling frequency must satisfy:

fs > 2fmax

This is the Nyquist condition. In real instruments, an analogue anti-aliasing filter is required because noise and interference above the intended bandwidth can otherwise fold into the measured band.

Every dataset should include metadata: sensor identity, calibration information, units, sampling rate, timestamp convention, gain, filter settings, and missing-value handling. A numerical array without these details is not an adequately defined engineering measurement.

5.2 Descriptive statistics

For N measurements x₁, …, x_N, the sample mean is:

x̄ = (1/N)∑i=1Nxi

The sample variance is:

s² = [1/(N−1)]∑i=1N(xi − x̄)²

and the sample standard deviation is s = √s². The N−1 denominator provides an unbiased estimator of population variance under independent, identically distributed observations with finite variance. It is not a correction for every type of dependence or measurement bias.

The root-mean-square value is:

xRMS = √[(1/N)∑i=1Nxi²]

For a zero-mean sinusoid, RMS equals peak magnitude divided by √2. For a signal with nonzero mean, RMS includes both the mean and fluctuating components.

A useful decomposition is:

xRMS² = x̄² + xac,RMS²

where x_ac,RMS is the RMS of x_i − x̄. This identity is exact for the finite dataset when consistent averaging is used.

5.3 Filtering and spectral analysis

A moving average of window length M is:

y[n] = (1/M)∑k=0M−1x[n−k]

It can reduce random noise but also smooth edges, introduce delay, and attenuate genuine high-frequency information. Filter choice must follow the engineering question rather than visual appearance alone.

The discrete Fourier transform (DFT) of N samples is:

X[k] = ∑n=0N−1x[n]e−j2πkn/N

where k = 0, …, N−1. For sampling frequency f_s, the frequency associated with bin k is kf_s/N for the unshifted DFT. Spectral leakage occurs when the sampled record does not contain an integer number of cycles of a periodic component. Window functions can reduce leakage but alter amplitude and noise bandwidth; amplitude correction may therefore be necessary.

6. Visualization and Communication

A plot is an analytical instrument, not decoration. Time-domain plots should label time units, signal units, sampling conditions, and relevant limits. Frequency plots should state whether magnitude is linear or decibel scale and identify the frequency resolution.

For a voltage gain G, decibels are:

GdB = 20 log₁₀(|Vout/Vin|)

For a power ratio P_out/P_in:

GdB = 10 log₁₀(Pout/Pin)

The different coefficients arise because power is proportional to squared voltage or current under consistent impedance conditions. Never use the power formula for an arbitrary voltage ratio without considering impedance.

Good engineering visualization includes:

  • meaningful axis labels and SI units;
  • a caption stating what was calculated;
  • legends only when multiple series require them;
  • uncertainty bands or error bars when measurement uncertainty matters;
  • avoidance of misleading truncated axes;
  • reproducible plot-generation code;
  • explicit treatment of missing, clipped, or saturated samples.

7. Floating-Point Error and Numerical Stability

7.1 Sources of error

Three broad error categories should be distinguished:

  1. Measurement error: uncertainty or bias in physical observations.
  2. Truncation error: error caused by replacing an exact mathematical process with an approximation, such as a finite difference or trapezoidal rule.
  3. Round-off error: error caused by finite-precision arithmetic.

Floating-point arithmetic is generally not associative. In exact arithmetic, (a + b) + c = a + (b + c); in floating point, the results may differ. Therefore, comparing computed values with == is unsafe when an approximate equality is intended.

A practical test uses absolute and relative tolerances:

|x − y| ≤ atol + rtol|y|

Here atol and rtol are application-dependent absolute and relative tolerances. They must be selected from physical accuracy requirements, not copied blindly from a software default.

7.2 Cancellation, overflow, and underflow

Catastrophic cancellation occurs when nearly equal numbers are subtracted, causing significant digits to be lost. For example, calculating a small difference between two large sensor readings can be much less accurate than measuring the difference directly.

Overflow occurs when a result exceeds the representable range; underflow occurs when it becomes too small for normal representation. Scaling variables, using stable rearrangements, and checking ranges can reduce these risks.

7.3 Stability and conditioning

Conditioning describes sensitivity of the mathematical problem; stability describes how much an algorithm amplifies computational errors. A well-conditioned problem can be ruined by an unstable algorithm, while a stable algorithm cannot remove sensitivity inherent in an ill-conditioned problem.

Best practices include:

  • scale variables so their magnitudes are comparable;
  • use library factorizations rather than explicit inverses;
  • avoid normal equations when QR or SVD is available;
  • inspect condition estimates for sensitive systems;
  • verify residuals and physical constraints;
  • repeat calculations at higher precision or with finer resolution when practical;
  • test limiting cases and conservation laws;
  • retain units and significant figures throughout the workflow.

8. Integrated Case Study: Monitoring a DC Motor Test

The following is an illustrative case, not a report of a particular commercial installation.

Background and problem

A laboratory team measures motor voltage, current, and shaft speed during a 10 s acceleration test. The objective is to estimate electrical energy consumed, identify the steady-state speed, and detect whether the current sensor contains a bias.

Analysis and approach

The data logger samples all channels at 1 kHz. Engineers first verify synchronized timestamps, units, sensor calibration, and missing samples. A plot of current against time reveals a constant offset of approximately 0.08 A during a no-load zero-current calibration interval. The corrected current is:

icorrected(t) = imeasured(t) − 0.08 A

Instantaneous electrical power is calculated as:

p(t) = v(t)icorrected(t)

The trapezoidal rule estimates energy:

E ≈ ∫010p(t)dt

The team computes the mean and standard deviation of speed over the final one-second window rather than over the entire acceleration period. This prevents transient behavior from being misreported as operating variability. A moving average is used only for display; the energy calculation uses the unfiltered, calibrated measurements because filtering would alter the signal and potentially bias the integral.

A root-finding routine identifies the first time at which speed reaches 95% of its final-window mean. Bisection is applied to a continuous interpolation of the sampled speed data over a bracket known to contain the crossing.

Outcome

The resulting report includes energy in joules, final speed in revolutions per minute, the standard deviation of the final-window speed, the calibration correction, sampling rate, integration method, and uncertainty limitations. The report also states that electrical energy is not equivalent to mechanical output energy because motor losses and mechanical load are not separately measured.

Learning points

The case demonstrates that numerical programming is a chain of reasoning rather than a collection of isolated commands. Calibration precedes statistics, units precede interpretation, and algorithm choice depends on the physical question. It also shows why a visually smooth curve should not automatically replace raw data in quantitative calculations.

9. Reproducible and Responsible Numerical Practice

A professional project should separate data acquisition, cleaning, computation, visualization, and reporting. Store raw data as read-only files and save processed data with a clear transformation record. Use version control for code, record package versions, and generate figures from scripts rather than editing values manually.

Numerical results may influence equipment ratings, protective settings, medical instrumentation, or energy decisions. Engineers must therefore document assumptions, validate code with known solutions, test edge cases, and report uncertainty honestly. AI-assisted programming can suggest code or identify possible errors, but generated code must be reviewed, tested, and traced to authoritative mathematical and engineering requirements. Confidential measurement data must not be uploaded to services without authorization.

Challenges and Limitations

Numerical methods approximate physical reality. A highly precise calculation based on an inaccurate sensor remains inaccurate. Sampling can hide fast transients; filtering can remove meaningful events; statistical summaries can conceal nonstationarity; and a small residual can coexist with a badly conditioned system. Numerical libraries also differ in defaults, handling of missing values, complex conventions, and tolerance policies. Reproducibility requires recording these details.

Best-Practice Checklist

  1. Define variables, units, sign conventions, and sampling assumptions.
  2. Validate input dimensions, ranges, timestamps, and missing values.
  3. Prefer stable library solvers to explicit matrix inversion.
  4. Check residuals, limiting cases, and physical conservation laws.
  5. Distinguish measurement, truncation, and round-off error.
  6. Use appropriate tolerances rather than exact equality for floating-point results.
  7. Plot raw and processed data when filtering affects interpretation.
  8. Record calibration, code version, parameters, and algorithm choices.
  9. Report uncertainty and limitations with the result.
  10. Protect data integrity, confidentiality, and safety-critical decisions.

Chapter Summary

Numerical and data programming provides the computational foundation for electrical engineering analysis. Vectors and matrices represent measurements and circuit models; reliable linear-system solvers support nodal, state-space, and estimation problems. Complex arithmetic enables phasor and impedance calculations. Numerical integration estimates accumulated quantities such as energy, while root-finding algorithms locate operating points and threshold crossings. Sampling theory, filtering, DFT analysis, and statistical summaries support signal and measurement processing.

Correct numerical practice requires more than obtaining a number. Engineers must evaluate units, assumptions, sampling conditions, conditioning, algorithmic stability, floating-point error, uncertainty, and physical plausibility. Reproducible code, documented metadata, residual checks, and honest visualization convert numerical output into defensible engineering evidence.

Key Terms

Floating-point arithmetic: Finite-precision representation and computation of numerical quantities.

Vector: An ordered one-dimensional collection of values.

Matrix: A rectangular array of numerical elements.

Condition number: A measure of the sensitivity of a problem to perturbations in its input.

Residual: The difference b − Ax after an approximate solution is substituted into a linear system.

Complex number: A quantity of the form a + jb, containing real and imaginary parts.

Phasor: A complex representation of a sinusoidal steady-state quantity.

Impedance: The complex ratio of voltage phasor to current phasor.

Numerical integration: Approximation of a definite integral using finite numerical data or evaluations.

Trapezoidal rule: An integration method that approximates intervals by trapezoids.

Root finding: The numerical determination of a value satisfying f(x) = 0.

Sampling frequency: The number of samples acquired per unit time.

Aliasing: Misrepresentation of frequency caused by insufficient sampling.

Discrete Fourier transform: A finite-sum transformation from sampled data to discrete frequency components.

RMS value: The square root of the mean of squared signal values.

Variance: A measure of squared dispersion about a mean.

Round-off error: Error caused by finite-precision numerical representation.

Truncation error: Error caused by replacing an exact process with an approximation.

Numerical stability: The degree to which an algorithm avoids amplifying computational errors.

Conditioning: The inherent sensitivity of a mathematical problem to input perturbations.

Review Questions

Short Answer Questions

  1. Why is solve(A, b) generally preferred to explicitly calculating A⁻¹b?
  2. State the impedance of an ideal inductor and define each symbol.
  3. What physical problem can aliasing cause in sampled measurements?
  4. Distinguish round-off error from truncation error.
  5. Why should measurement units and calibration metadata be retained with numerical data?

Long Answer Questions

  1. Explain how complex phasors simplify sinusoidal steady-state circuit analysis. Include the relationship between RMS phasors and time-domain sinusoids.
  2. Derive or explain the composite trapezoidal rule and discuss conditions that affect its accuracy.
  3. Compare bisection and Newton’s methods for root finding, including assumptions, convergence behavior, and failure modes.
  4. Explain the difference between conditioning and numerical stability, and describe how each affects the solution of a linear system.
  5. Design a reproducible workflow for processing voltage and current measurements from an electrical experiment.

Analytical / Application Questions

  1. A sampled power record contains voltage in volts and current in amperes at nonuniform timestamps. Describe how you would calculate energy, validate the result, and report its units and limitations.
  2. A matrix-based circuit model produces a small residual but an extremely large condition number. Explain what this indicates and recommend additional checks before accepting the node voltages.
  3. A motor-speed signal is noisy and includes a startup transient. Propose a method for estimating steady-state speed and its variability without allowing the transient or display filtering to distort the engineering conclusion.

Multiple Choice Questions with Answers

Question 1

Which operation is normally preferred for solving A x = b in numerical engineering software?

A. Compute the explicit inverse of A and multiply by b

B. Use a numerical linear-system solver such as LU-based solve

C. Add all entries of A and divide b by the result

D. Use elementwise multiplication between A and b

Correct Answer: B. Use a numerical linear-system solver such as LU-based solve

Question 2

For a sinusoidal steady-state inductor, which impedance is correct?

A. R + jωL

B. 1/(jωL)

C. jωL

D. 1/(jωC)

Correct Answer: C. jωL

Question 3

A continuous function has f(a)f(b) < 0. Which root-finding method directly uses this bracket condition?

A. Bisection

B. Unbracketed Newton iteration only

C. Moving average

D. Discrete Fourier transform

Correct Answer: A. Bisection

Question 4

If a signal has a highest relevant frequency of 4 kHz, which sampling frequency satisfies the ideal Nyquist inequality?

A. 4 kHz

B. 6 kHz

C. 8 kHz

D. 10 kHz

Correct Answer: D. 10 kHz

Question 5

What does a large condition number generally indicate?

A. The problem is insensitive to input errors

B. Small input perturbations may cause large solution changes

C. The matrix necessarily has a zero determinant

D. The residual must be exactly zero

Correct Answer: B. Small input perturbations may cause large solution changes

Question 6

For N independent measurements, which denominator is used in the usual sample variance?

A. N−1

B. N+1

C. N²

D. 2N

Correct Answer: A. N−1

Question 7

Which expression gives the RMS value of samples x_i?

A. (1/N)∑x_i

B. √[(1/N)∑x_i²]

C. ∑|x_i|/N²

D. √[∑x_i/N]

Correct Answer: B. √[(1/N)∑x_i²]

Question 8

Why can explicit formation of AᵀA be undesirable in least-squares computation?

A. It always changes real data into complex data

B. It makes the number of measurements zero

C. It can worsen the conditioning of the problem

D. It prevents any residual from being calculated

Correct Answer: C. It can worsen the conditioning of the problem

Question 9

For a voltage ratio measured under consistent impedance conditions, which decibel expression is appropriate?

A. 10 log₁₀(|V_out/V_in|)

B. 20 log₁₀(|V_out/V_in|)

C. log₂(|V_out/V_in|)

D. (|V_out/V_in|)²

Correct Answer: B. 20 log₁₀(|V_out/V_in|)

Question 10

Which quantity should be checked to assess how well an approximate x satisfies A x = b?

A. The residual b − Ax

B. The number of plot colors

C. The sampling frequency alone

D. The determinant of b

Correct Answer: A. The residual b − Ax

Programming Fundamentals for Engineers
Chapter 13 of 15