Discretization: How to Program a PID in C++

The discretization is the process of converting continuous equations into discrete steps so they can be executed on a microcontroller or computer.

In theoretical articles, we’ve casually talked about “integrating the error” and “calculating the derivative.” The problem is that the microcontrollers we commonly use, like those in an Arduino or ESP32, have no idea what differential calculus is. They only know how to add, subtract, multiply, and divide numbers in a sequential loop that repeats over and over.

Today, we’re coming down from the clouds of Laplace theory and getting our hands dirty with silicon. We’re going to see how to convert that math into clean, efficient C++ code.

From Continuous Model to Digital Model

To refresh your memory, the classic PID controller equation in the continuous time domain is:

In this equation, time t flows smoothly and constantly. But in our microcontroller, time advances “in jumps.” We read the sensor, calculate, write to the actuator, and start over. Each of these cycles is called the sampling time or Δt (dt).

Our mission is to replace the continuous function e(t) with an instantaneous reading at a specific moment e[k], where k is the current loop cycle, and k−1 is the previous cycle.

Proportional Action (P)

This is the easy part. Proportional action has no memory; it only depends on the present moment. Therefore, its translation into code is straightforward. We simply multiply the current error by the gain:

Integral Action (I): Summing Rectangles

Here’s where things get interesting. The integral of a signal is, by geometric definition, the area under the curve that the signal forms.

Since our microcontroller can’t see the perfect curve, it only sees a new data point every Δt milliseconds. To calculate the area, we use a numerical technique called Euler Integration. It involves assuming that during time Δt, the error remained constant, forming a small rectangle with base Δt and height e[k].

The total accumulated area is simply the area we already had stored plus the area of this new rectangle:

Derivative Action (D): Calculating Slopes

The derivative represents the rate of change or the slope of the curve at a given instant.

In our discrete-step world, if we want to know how fast the error is changing, we simply subtract the current error minus the error from the previous cycle, and divide it by the time that has passed between the two measurements (again, our beloved Δt):

For this to work in code, at the end of each cycle, we must remember to do previous_error = current_error, saving it for the future.

Implementation in C++ for Arduino/ESP32

Now, let’s put this into real code. To calculate the elapsed time (Δt), in the Arduino ecosystem we usually use the millis() function (or micros() for very fast systems).

A tempting solution would be to use a delay(10). In a program that shares the processor with communications, interfaces, or other tasks, it’s best to avoid delay() inside the control loop, because it blocks its execution and makes it difficult to control the actual period.

Let’s create a basic structure using the “timestamp evaluation” paradigm:

// PID tuning global variables
float Kp = 2.0;
float Ki = 5.0;
float Kd = 1.0;

// PID internal state variables (Memory)
float accumulated_integral = 0.0;
float previous_error = 0.0;
unsigned long previous_time = 0;

// System variables
float setpoint = 100.0;

void setup() {
    // Hardware initialization
    previous_time = millis();
}

void loop() {
    unsigned long current_time = millis();
    unsigned long dt_ms = current_time - previous_time;

    // Only execute the calculation if a minimum time has passed
    // For example, run the PID every 10ms (100 Hz)
    if (dt_ms >= 10) {
        
        // 1. Convert dt to seconds (essential for Ki and Kd constants to have physical meaning)
        float dt = (float)dt_ms / 1000.0; 
        
        // 2. Measure the process variable
        float pv = readSensor();
        
        // 3. Calculate the error
        float error = setpoint - pv;
        
        // --- PID CALCULATION ---
        
        // Proportional
        float P = Kp * error;
        
        // Integral (Euler integration)
        accumulated_integral += error * dt;
        float I = Ki * accumulated_integral;
        
        // Derivative (Finite difference)
        float derivative = (error - previous_error) / dt;
        float D = Kd * derivative;
        
        // Total output
        float output = P + I + D;
        
        // --- END OF CALCULATION ---
        
        // 4. Actuate on the plant
        writeActuator(output);
        
        // 5. Save state for the next cycle
        previous_error = error;
        previous_time = current_time;
    }
    
    // Other non-blocking tasks can go here
    // read_serial_port();
    // update_display();
}

Critical Notes on the Implementation

1. The conversion to seconds: Notice that we divide dt_ms by 1000.0. This is important. If you don’t do this and multiply directly by milliseconds, your integral term will become enormous in fractions of a second, and your derivative will be tiny, forcing you to use weird, unstable Ki and Kd values.
2. Determinism: The approach with millis() in the loop() is fine for basic robotics. However, if your microcontroller is doing many things (like handling a web server with WiFi on an ESP32), the dt_ms time could fluctuate between 10 ms and 15 ms. If you need to reduce jitter, use a hardware timer or a periodic task. Keep the ISR brief and leave the calculation and peripheral access in a task with controlled priority and execution time.

The Elephant in the Room

The code we just wrote reproduces the basic mathematical structure, but it’s not ready for a real plant yet. It’s missing output limits, anti-windup, derivative filtering, and explicit handling of the first cycle to avoid a derivative spike.

The real world also introduces several problems. What happens if we calculate that the motor needs “2000 Volts” but our battery only provides 12V? The integral will start accumulating a gigantic error that it can never erase. And what if our sensor has electrical noise? Then the subtraction e[k]−e[k−1] in the derivative will generate wild calculation spikes that will make our motor screech and vibrate like a broken washing machine.

These are the real problems of digital implementation, and they are exactly what we are going to solve in the last two articles of the course. Starting with the dreaded integral “Windup.”