Simulating 2nd Order System with ESP32

Question

How to simulate 2nd order system with ESP32. the system has 1 analog input and 1 analog output.

Answer

To simulate a second-order system using an ESP32 with one analog input and one analog output, follow these steps:

Step-by-Step Explanation

Understand the System: The goal is to simulate a continuous-time second-order system discretized in time. This involves reading an analog input, processing it through a discrete model of the system, and sending the result to an analog output.
Discretize the System: Convert the continuous-time second-order differential equation into a discrete form using methods like the bilinear transform or forward/backward Euler. Here, we use the discretized state-space representation derived from the transfer function ωn2s2+2ζωns+ωn2s2+2ζωns+ωn2ωn2, resulting in the difference equation: yk+1=a1yk+a2yk−1+bukyk+1=a1yk+a2yk−1+buk where a1=2−2ζωnT−ωn2T2a1=2−2ζωnT−ωn2T2, a2=−1+2ζωnTa2=−1+2ζωnT, and b=ωn2T2b=ωn2T2.
Setup the Hardware: Connect an analog sensor to the ESP32’s ADC pin (e.g., GPIO36) and connect a DAC output (e.g., GPIO25) to drive the system’s response.
Read Analog Input: Use the ESP32’s ADC to read the input voltage, scaling it to a normalized range [0,1].
Compute Output: Apply the discrete model using the current input and previous outputs to compute the new output value.
Scale and Clamp: Ensure the computed output stays within the DAC’s voltage range by scaling and clamping.
Output the Result: Convert the computed output into a DAC value and write it to the analog pin.
Update States: Keep track of the previous two outputs for the next iteration.

Example Code

#include <Arduino.h>
#include "driver/dac.h"

#define ADC_PIN GPIO_NUM_36   // Analog input pin
#define DAC_CHANNEL DAC_CHANNEL_0  // DAC channel for output

void setup() {
    dac_config_t config;
    dac_default_config(&config);
    config.dac_num = DAC_CHANNEL_0;
    config.resolution = DAC_RESOLUTION_8_BIT; // Use 8-bit resolution
    dac_config(&config);

    // Enable the selected DAC channel
   dac_enable(DAC_CHANNEL_0, true);
}

void loop() {
    // Read ADC value
    int adc_val = analogRead(ADC_PIN);
    
    // Parameters for the second-order system
    const float zeta = 0.5;      // Damping ratio
    const float wn = 10.0;       // Natural frequency in rad/s
    const float T = 0.001;       // Sampling period (1ms)
    
    static float y_prev = 0.0;
    static float y_pprev = 0.0;

    // Compute coefficients
    float a1 = 2 - 2*zeta*wn*T - wn*wn*T*T;
    float a2 = -1 + 2*zeta*wn*T;
    float b = wn*wn*T*T;

    // Normalize ADC input to [0,1]
    float u = (float)adc_val / 4095.0; 

    // Compute current output
    float y_current = a1 * y_prev + a2 * y_pprev + b * u;

    // Scale y_current to DAC voltage range [0,3.3V] with 8-bit resolution
    if (y_current < 0) y_current = 0;
    else if (y_current > 3.3) y_current = 3.3;
    
    int dac_val = static_cast<int>((y_current / 3.3) * 255);

    // Write to DAC
    dac_write(DAC_CHANNEL_0, dac_val);

    // Update states for next iteration
    y_pprev = y_prev;
    y_prev = y_current;

    // Wait for the sampling period
    delay(1);
}

Explanation of Code

Initialization: The setup function configures the DAC with 8-bit resolution on channel 0.
Loop Execution:
- Reads the analog input from ADC pin 36.
- Defines system parameters: damping ratio (ζ), natural frequency (ωn), and sampling period (T).
- Computes the coefficients a1, a2, and b based on the discretized model.
- Normalizes the ADC value to [0,1].
- Calculates the current output using the state equation.
- Scales and clamps the output to fit within the DAC’s voltage range.
- Writes the computed value to the DAC.
- Updates previous states for the next iteration.
Timing: Uses delay(1) to maintain a 1ms sampling period, ensuring consistent discretization.

This approach allows real-time simulation of a second-order system on an ESP32, enabling experimentation with different damping ratios and natural frequencies while observing the system’s response through analog output.

Thinking Process

Okay, so I’m trying to figure out how to simulate a second-order system using an ESP32. The system has one analog input and one analog output. Hmm, I remember that second-order systems are often described by differential equations, maybe something like a mass-spring-damper system or an RLC circuit. But since this is about simulation with microcontrollers, I think the focus will be more on implementing the mathematical model of such a system.

First, I need to recall what a second-order system typically looks like in terms of its transfer function. It’s usually expressed as G(s) = (ωn²)/(s² + 2ζωn s + ωn²), where ζ is the damping ratio and ωn is the natural frequency. So, maybe I need to implement this transfer function on the ESP32.

But wait, how do you translate that into something an ESP32 can handle? The ESP32 can take analog inputs via ADC pins and output analog signals through DAC pins. So, the input will be read from an analog pin, processed in some way, and then sent out as an analog signal.

I remember that for digital implementations of continuous systems, we often use numerical integration methods like Euler or Runge-Kutta. Maybe I can discretize the second-order differential equation so that it can be computed step by step on the microcontroller.

Let me think about the differential equations. A second-order system typically has two derivatives. So if I have a system described by d²y/dt² + 2ζωn dy/dt + ωn² y = ωn² u(t), where u(t) is the input, and y(t) is the output.

To convert this into discrete-time equations, I can use Euler’s method. That means approximating derivatives with finite differences. So, the second derivative d²y/dt² ≈ (y_new – 2*y_old + y_older)/T², where T is the sample time.

Similarly, the first derivative dy/dt ≈ (y_new – y_old)/T.

Wait, but when implementing this, I’ll need to keep track of previous states. So maybe I should set up variables for y(t), y(t-1), and y(t-2). Similarly for the input u(t).

So, rearranging the equation:

y_new = (2ζωn T²) * (y_old – y_older) + ωn² T² u_current

Wait, no, maybe I need to think about it differently. Let me write out the discretized version step by step.

The original differential equation is:

d²y/dt² = -2ζωn dy/dt – ωn² y + ωn² u(t)

Using Euler’s method for time discretization, let T be the sample period.

Let’s denote y_k as y(k*T), and similarly for other variables.

Then,

(y_{k+1} – 2y_k + y_{k-1}) / T² = -2ζωn (y_k – y_{k-1})/T – ωn² y_k + ωn² u_k

Multiplying both sides by T²:

y_{k+1} – 2y_k + y_{k-1} = -2ζωn T (y_k – y_{k-1}) – ωn² T² y_k + ωn² T² u_k

Then, solving for y_{k+1}:

y_{k+1} = 2y_k – y_{k-1} – 2ζωn T (y_k – y_{k-1}) – ωn² T² y_k + ωn² T² u_k

Simplify terms:

Let me expand that:

y_{k+1} = [2y_k] + [- y_{k-1}] + [-2ζωn T y_k + 2ζωn T y_{k-1}] + [-ωn² T² y_k] + [ωn² T² u_k]

Group like terms:

Terms with y_k: 2y_k -2ζωn T y_k – ωn² T² y_k

Terms with y_{k-1}: -y_{k-1} + 2ζωn T y_{k-1}

Plus the term with u_k: ωn² T² u_k

So,

y_{k+1} = (2 – 2ζωn T – ωn² T²) y_k + (-1 + 2ζωn T) y_{k-1} + ωn² T² u_k

That’s the recurrence relation I’ll need to implement.

So, on the ESP32, at each sample time T, I’ll read the current input u_k from the analog pin, then compute y_{k+1} using this equation. Then, output y_{k+1} through the DAC pin.

I’ll also need to store the previous two outputs (y_k and y_{k-1}) as state variables.

But how do I handle the sample time T? Since the ESP32 can run loops at a certain frequency, I can set up a timer interrupt or use a loop with a delay. But delays in loops are generally bad for responsiveness, so maybe using an interrupt-based approach is better. Alternatively, using the built-in functions like delay() but making sure that the sampling period is consistent.

Let me outline the steps:

Set up the ADC and DAC pins. The analog input will be read from a pin, say GPIO36 (ADC0), and the output will go through a DAC pin, maybe GPIO25 (DAC0).
Initialize variables: y_prev_prev = 0, y_prev = 0.
Read the current input u from ADC.
Compute y_current using the recurrence relation:y_current = a1 * y_prev + a2 * y_prev_prev + b * u

Where: a1 = (2 – 2ζωn T – ωn² T²) a2 = (-1 + 2ζωn T) b = ωn² T²

Update the DAC output with y_current.
Shift the state variables: y_prev_prev becomes y_prev, and y_prev becomes y_current.
Repeat at intervals of T seconds.

Wait, but in code terms, this would be a loop that runs every T seconds, reads the input, computes the next output based on previous states and current input, then outputs it.

Potential issues:

Choosing appropriate sampling time T: It should be small enough to capture the dynamics accurately without causing numerical instability. The Nyquist rate might give some guidance.
Scaling of ADC and DAC values: ESP32’s ADC reads 12-bit values (0-4095), DAC outputs also 8 or 12 bits, depending on configuration. So I’ll need to scale the computed y_current appropriately to fit into the DAC output range.

For example, if y_current is a float between -1 and 1, but DAC expects 0-255 (for 8-bit), then after computing y_current, we can map it to the DAC range. Or maybe just use the entire ADC input range as 0-3.3V, so scaling accordingly.

Also, I might need to filter the ADC input to avoid noise, perhaps using a moving average or low-pass filter, but that could complicate things. For simulation purposes, maybe it’s okay without filtering initially.

Another thought: Should I use floating-point arithmetic? On an ESP32, which is a 32-bit processor, yes, floats should be manageable, but I need to make sure the computations are precise enough.

Let me think about the code structure:

Include necessary headers for ADC and DAC.

Initialize ADC with appropriate settings (e.g., input pin, sample rate). Similarly, initialize DAC output.

Set up variables:

float y_prev = 0.0; float y_pprev = 0.0;

Parameters ζ and ωn are given or chosen based on the desired system behavior. T is the sampling period.

Wait, but in code, how do I manage timing? Maybe set up a timer interrupt every T seconds that triggers the computation and output. Alternatively, use a loop with delay(T). Using delay() could cause issues if other tasks need to run, but for a simple simulation, it might suffice.

But using an interrupt would be better for real-time performance, ensuring that computations happen exactly at each sampling interval without blocking other operations.

So perhaps setting up Timer 0 or another timer to trigger every T milliseconds.

In the interrupt service routine (ISR), read ADC, compute y_current, write DAC, update states.

Wait, but the ADC reading might take some time. So I need to make sure that the ISR doesn’t block for too long. Alternatively, set up the ADC in a way that allows non-blocking reads or use DMA for ADC transfers. But this could get complicated.

Alternatively, read the ADC value outside the ISR and store it so that the ISR can quickly access it. Hmm, but timing might be tricky.

Wait, perhaps I should use a loop with delay(T) in the main function. Since the ESP32 is single-threaded (for user code), using delays will block other operations, but if this is the only task, maybe it’s acceptable.

So code outline:

#include <ADC.h> #include <DAC.h>

// Initialize ADC and DAC ADC adc; DAC dac;

void setup() { // Initialize ADC on pin 36 with default settings. adc.config(); // Initialize DAC on pin 25 (or another DAC pin), set output mode to voltage, and resolution. dac.enable(DAC_CHANNEL_0); }

void loop() { // Read ADC value int adc_val = adc.read(ADC_PIN_36); // Assuming pin 36 is used.

// Compute y_current based on previous states and current input
float u = (float)adc_val / 4095.0;  // Scale to [0,1]

// Parameters: ζ, ωn, T
const float zeta = 0.5;
const float wn = 2.0;
const float T = 0.01;  // Sampling period in seconds

static float y_prev = 0.0;
static float y_pprev = 0.0;

float a1 = (2 - 2*zeta*wn*T - wn*wn*T*T);
float a2 = (-1 + 2*zeta*wn*T);
float b = wn*wn*T*T;

float y_current = a1 * y_prev + a2 * y_pprev + b * u;

// Scale y_current to DAC range
int dac_val = (int)(y_current * 255);  // Assuming 8-bit DAC, scale from [-1,1] or [0,1]? Wait, need to adjust.

// If y_current is in a certain range, map it appropriately.
// For example, if y ranges from -1 to 1:
dac_val = (int)((y_current + 1.0) * 127.5);  // So that -1 maps to 0 and 1 maps to 255.

// Ensure dac_val is within 0-255
if (dac_val < 0) dac_val = 0;
else if (dac_val > 255) dac_val = 255;

// Write DAC value
dac.write(DAC_CHANNEL_0, dac_val);

// Update states
y_pprev = y_prev;
y_prev = y_current;

// Wait for next sampling interval
delay(1);  // Assuming T is 1ms? Or adjust based on desired T.

Wait, but in the code above, the computation of a1, a2, and b is done every loop. But these are constants once ζ, ωn, and T are set. So maybe precompute them outside the loop to save processing time.

Also, in the example, I assumed y_current could be between -1 and 1, but depending on the system’s behavior, it might not reach those limits. Alternatively, scale the input u appropriately if needed.

Wait, another thing: The initial conditions are set to y_prev = 0 and y_pprev = 0. Depending on the simulation, maybe they should be initialized differently.

Testing this code would involve applying a step input (e.g., connecting VCC to ADC pin after some time) and observing the output waveform. If it behaves like a second-order system with the chosen ζ and ωn, then it’s working.

Potential improvements:

Use higher-resolution DAC by setting it to 12-bit mode if available.
Implement anti-aliasing filters on the ADC input to reduce noise.
Add more sophisticated state storage, maybe using arrays for a better state representation.
Consider using floating-point math carefully to avoid overflow or underflow.

Wait, but in my code outline above, I have:

float u = (float)adc_val / 4095.0;

Assuming ADC is 12-bit. Then, scaling to [0,1]. But depending on how the system is set up, maybe we want u to be in a different range.

Similarly, when computing y_current, if it’s supposed to be in [-1,1], but due to the coefficients, it might go beyond that. So scaling needs to be handled carefully.

Another consideration: The choice of T affects both the discretization and the computation rate. If T is too large, the simulation may not capture high-frequency dynamics accurately and could introduce errors or even instability if the sampling rate is too low relative to the system’s bandwidth.

So maybe choosing a smaller T, like 1ms (T=0.001s), which would correspond to a sample rate of 1kHz. That should be sufficient for many second-order systems unless they have very high natural frequencies.

Wait, in my earlier code example, I set delay(1) which is 1ms, so T=0.001s.

Another thing: The computation inside the loop needs to finish before the next sampling interval; otherwise, it could cause delays and inconsistent timing. So make sure that the computations are quick enough for the chosen T.

In terms of testing, I might want to log the output values over time and plot them to see if they match the expected second-order response (like overdamped, underdamped, critically damped).

Alternatively, use an oscilloscope to observe the DAC output when a step input is applied via the ADC.

Potential issues in code:

Incorrect scaling between ADC and DAC. For example, if the ADC reads 0-4095 for 0-3.3V, and the DAC outputs 0-255 for 0-3.3V (8-bit), then scaling should be consistent. So perhaps u = (adc_val / 4095.0) * 3.3; but if the DAC expects the same voltage range, then y_current is in volts, so it needs to be scaled accordingly.

Wait, maybe I should represent all variables in terms of voltage. Let me think:

ADC input: 0-3.3V mapped to 0-4095.

So u_volts = (adc_val / 4095) * 3.3;

DAC output: 0-255 maps to 0-3.3V (assuming 8-bit, which is the case for ESP32’s DAC).

Wait, no, ESP32’s DAC can be set in different resolutions. By default, it might be 8 bits.

So if I use an 8-bit DAC:

dac.write(DAC_CHANNEL_0, value); where value is 0-255.

Then, 0 corresponds to 0V, 255 corresponds to 3.3V.

Therefore, to set a voltage V_out, the DAC value should be (V_out / 3.3) * 255.

So in code:

float y_volts = …; // computed output in volts int dac_val = (int)(y_volts / 3.3 * 255);

But if y_volts can go beyond the 0-3.3V range, we need to clamp it.

So perhaps after computing y_current (which is in volts), scale and clamp:

if (y_volts < 0) y_volts = 0; else if (y_volts > 3.3) y_volts = 3.3;

Then compute dac_val as above.

Alternatively, design the system so that y remains within this range.

Wait, but in my previous code example, I didn’t consider voltage units; instead, I treated u as a normalized value between 0 and 1. That could lead to confusion because then y_current is also scaled that way. Maybe it’s better to keep everything in volts for clarity.

So let me adjust:

u_volts = (adc_val / 4095) * 3.3;

Then compute y_volts based on the differential equation, but ensuring that y_volts remains within [0,3.3] or maybe [-3.3,3.3], depending on the system.

Wait, in many second-order systems, especially in control applications, inputs and outputs are often zero-centered, so maybe it’s better to have u range from -1 to 1 (or -3.3V to +3.3V), but ADC can’t read negative voltages directly unless using an external op-amp or something.

Therefore, perhaps the input is DC shifted to be within 0-3.3V, and similarly for the output.

Alternatively, use a bipolar range by adding an offset. For example:

ADC reads 1.65V as zero, so subtracting 1.65V gives u_volts in -1.65V to +1.65V. Similarly for DAC output: add 1.65V and scale accordingly.

But this adds complexity. Maybe it’s better to stick with unipolar for simplicity unless the simulation requires a zero-centered input/output.

So, perhaps proceed with ADC reading as [0-3.3V] mapped to [0-4095], then scaled to [0-1]. Similarly, DAC output is scaled from [0-1] back to [0-3.3V].

Therefore, in code:

u = (float)adc_val / 4095.0; // u ranges [0,1]

Compute y_current, which should be a value between 0 and 1 as well.

But wait, depending on the system’s natural frequency and damping ratio, y_current could exceed these limits, leading to saturation. So maybe include some scaling factors or limiters.

Alternatively, perhaps define y_current in a certain range and then scale it appropriately for DAC output.

Wait, but the equation for y_{k+1} is linear, so if u is bounded between 0 and 1, depending on the coefficients a1, a2, b, y could potentially exceed these bounds. So adding some form of limiting might be necessary to prevent the DAC from being driven beyond its capabilities.

In summary, simulating a second-order system with ESP32 involves reading an analog input, applying the difference equation (discretized version) based on the chosen parameters ζ and ωn, then outputting the result via DAC. The key steps are:

Read ADC value.
Convert to voltage or normalized value.
Apply discrete-time model: y_current = a1y_prev + a2y_pprev + b*u
Scale and clamp y_current to fit DAC output range.
Write to DAC.
Update previous states.

Potential areas for error:

Incorrect scaling between ADC and DAC, leading to unexpected output ranges.
Choosing too large a T, causing discretization errors or instability.
Not handling state updates correctly (e.g., overwriting y_prev before using it).
Forgetting to static declare y_prev and y_pprev, causing them to reset each loop.

Testing would involve applying known inputs (like step, sine wave) and verifying the output’s frequency response or transient behavior matches expectations based on ζ and ωn.

Finally, for better performance, especially with faster systems, using interrupts or more precise timing methods than delay() might be necessary. However, for many applications, delay() is sufficient.