Click here to Skip to main content
15,897,518 members
Articles / Programming Languages / C#

Slaving Over a Hot Laptop,PID control loops and Sous Vide cooking

Rate me:
Please Sign up or sign in to vote.
4.81/5 (12 votes)
12 Feb 2012CPOL10 min read 40.1K   1K   25  
PID thermal control loop simulation
///////////////////////////////////////////////////////////////////////////////
//
//  Form1.cs
//
//  By Philip R. Braica (HoshiKata@aol.com, VeryMadSci@gmail.com)
//
//  Distributed under the The Code Project Open License (CPOL)
//  http://www.codeproject.com/info/cpol10.aspx
///////////////////////////////////////////////////////////////////////////////

// Using.
using System;
using System.Collections.Generic;
using System.ComponentModel;
using System.Data;
using System.Drawing;
using System.Linq;
using System.Text;
using System.Windows.Forms;
using System.Windows.Forms.DataVisualization.Charting;

namespace SousVide
{
    public partial class Form1 : Form
    {
        public Form1()
        {
            InitializeComponent();
        }

        public void RunSimulation(bool usePID)
        {
            Parameters p = new Parameters(this);
            PID pid = new PID(p.Goal, p.P, p.I, p.D);
            double water = p.WaterStart;
            double container = p.WaterStart;
            double sensor = p.WaterStart;
            double food = p.WaterStart;
            double dt = 1.0 / p.Hz;

            Series s = new Series("Sensor");
            Series c = new Series("Container");
            Series w = new Series("Water");
            Series f = new Series("Food");

            Series we = new Series("Water Error");
            Series se = new Series("Measured Error");
            Series fe = new Series("Food Error");
            
            s.ChartType = System.Windows.Forms.DataVisualization.Charting.SeriesChartType.Line;
            c.ChartType = System.Windows.Forms.DataVisualization.Charting.SeriesChartType.Line;
            w.ChartType = System.Windows.Forms.DataVisualization.Charting.SeriesChartType.Line;
            f.ChartType = System.Windows.Forms.DataVisualization.Charting.SeriesChartType.Line;
            we.ChartType = System.Windows.Forms.DataVisualization.Charting.SeriesChartType.Line;
            se.ChartType = System.Windows.Forms.DataVisualization.Charting.SeriesChartType.Line;
            fe.ChartType = System.Windows.Forms.DataVisualization.Charting.SeriesChartType.Line;

            s.Color = Color.Black;
            c.Color = Color.Blue;
            w.Color = Color.Green;
            f.Color = Color.Orange;

            se.Color = Color.Black;
            we.Color = Color.Green;
            fe.Color = Color.Orange;

            // To make the graph reasonable, limit the graph points by downsampling.
            // The number samples without downsampling is p.Duration / dt.
            // Suppose we wanted 5k points, and p.Duration/dt is 15k, then the 
            // graphLimit would need to be 3 dt. 

            double graphLimit = p.Duration / 4000;
            double gt = 0;
            double adjt = 0;
            double maxOvershoot = 0;
            List<double> lastPart = new List<double>();
            double percentHeat = 0;
            double timeFirstWithinPointOne = 0;
            bool firstOverPointOne = false;
            for (double t = 0; t < p.Duration; t += dt)
            {
                // Throttle to the adjustment rate, assumes user is running simulation 10x faster than 
                // update rate.
                if (t - adjt > (1 / p.ControlRate))
                {
                    percentHeat = usePID ? pid.Update(sensor, dt) : p.HeatPercent;
                }
                simulateTime(dt, percentHeat, ref water, ref container, ref sensor, ref food, p, pid, usePID);
                double ft = food < -100 ? -100 : food > 1000 ? 1000 : food;
                if (firstOverPointOne == false)
                {
                    if (ft > p.Goal - 0.1)
                    {
                        timeFirstWithinPointOne = ft;
                        firstOverPointOne = true;
                    }
                }
                if (((t - gt) > graphLimit) || (t == 0))
                {                    
                    double st = sensor < -100 ? -100 : sensor > 1000 ? 1000 : sensor;
                    double wt = water < -100 ? -100 : water > 1000 ? 1000 : water;
                    s.Points.AddXY(t, st);
                    c.Points.AddXY(t, container < -100 ? -100 : container > 1000 ? 1000 : container);
                    w.Points.AddXY(t, wt);
                    f.Points.AddXY(t, ft);
                    if (ft > p.Goal)
                    {
                        maxOvershoot = maxOvershoot > ft ? maxOvershoot : ft;
                    }
                    if (t > p.Duration/2)
                    {
                        lastPart.Add(ft);
                    }
                    se.Points.AddXY(t, st - p.Goal);
                    we.Points.AddXY(t, wt - p.Goal);
                    fe.Points.AddXY(t, ft - p.Goal);
                    gt = t;
                }
            }

            chart1.Series.Clear();
            chart1.ChartAreas[0].RecalculateAxesScale();
            chart1.Series.Add(w);
            chart1.Series.Add(c);
            chart1.Series.Add(s);
            chart1.Series.Add(f);
            chart1.Update();

            chart2.Series.Clear();
            chart2.ChartAreas[0].RecalculateAxesScale();
            chart2.Series.Add(we);
            chart2.Series.Add(se);
            chart2.Series.Add(fe);
            chart2.Update();

            double ave = 0;
            double stdev = 0;
            if (lastPart.Count > 0)
            {
                for (int i = 0 ; i < lastPart.Count; i++)
                {
                    ave += lastPart[i];
                }
                ave /= lastPart.Count;
                for (int i = 0; i < lastPart.Count; i++)
                {
                    double err = lastPart[i] - ave;
                    stdev += (err*err);
                }
                stdev = System.Math.Sqrt(stdev)/lastPart.Count;
            }
            AnalysisLabel.Text =
                "Max overshoot = " + maxOvershoot.ToString("#0.000") + "\r\n" +
                (firstOverPointOne ? 
                ("Achieved 0.1 degree below goal @ " + timeFirstWithinPointOne.ToString("#0.0") + "s\r\n") :
                ("Never reached within 0.1 degree of goal.\r\n")) +
                "Ave. food tmp. last half of cooking = " + ave.ToString("#0.00") + "\r\n" +
                "Std. dev. of last half = " + stdev.ToString("#0.0000");
        }

        /// <summary>
        /// Convert watts to temperature, for now using 1.5.
        /// Need to tune to actual transfer using steady state analysis.
        /// </summary>
        private double heaterSpecificHeatFactor = 1.5;

        /// <summary>
        /// Update simulation.
        /// </summary>
        /// <param name="dt"></param>
        /// <param name="water"></param>
        /// <param name="container"></param>
        /// <param name="sensor"></param>
        /// <param name="p"></param>
        /// <param name="pid"></param>
        /// <param name="usePID"></param>
        private void simulateTime(double dt, double percentHeat, ref double water, ref double container, ref double sensor, ref double food, Parameters p, PID pid, bool usePID)
        {
            // Note no heat capacity terms are used because they are in the time constant after 
            // converting watts into transferable heat.

            // Limit to the max / min.
            percentHeat = percentHeat < 0 ? 0 : percentHeat > 100 ? 100 : percentHeat;

            // If binary, heater is on/off in which case threshold.
            if (p.Binary) { percentHeat = percentHeat > 25 ? 100 : 0; }

            double heaterAdded = percentHeat * p.HeaterWatts * heaterSpecificHeatFactor;

            double effectiveWaterContainer = p.TWC * p.WaterRadius * p.WaterDepth / (p.WaterRadius + 2 * p.WaterDepth);
            double effectiveWaterAir = p.TWA * p.WaterDepth; // surface area is constant but thermal mass of the water isn't.

            // Heat added from heater to water.
            double heatAddedToWater = (heaterAdded - water) * (1-System.Math.Exp(-dt / p.TWH))/p.TWH;
    heatAddedToWater = heatAddedToWater < 0 ? 0 : heatAddedToWater;

    // Heat lost to the sensor.
    double heatLostToSensor = (water - sensor) * (1-System.Math.Exp(-dt / p.TWS))/p.TWS;

    // Heat lost to the food.
    double heatToFood = (water - food) * (1 - System.Math.Exp(-dt / p.TWF)) / p.TWF;

    // Heat from water to container.
    double heatLostWaterContainer = (water - container) * (1 - System.Math.Exp(-dt / effectiveWaterContainer))/effectiveWaterContainer;

    // Heat water to air.
    double heatLostWaterAir = (water - p.Air) * (1 - System.Math.Exp(-dt / effectiveWaterAir))/effectiveWaterAir;

    // Heat Container to air.
    double heatLostContainerAir = (container - p.Air) * (1 - System.Math.Exp(-dt / p.TCA))/p.TCA;

    // Update.
    food = food + heatToFood;
    water = water + heatAddedToWater - heatLostWaterContainer - heatLostWaterAir - heatLostToSensor - heatToFood;
    sensor = sensor + heatLostToSensor;
    container = container + heatLostWaterContainer - heatLostContainerAir;
        }





        /// <summary>
        /// Grab GUI values into a bag.
        /// </summary>
        public class Parameters
        {
            /// <summary>
            /// Constructor.
            /// </summary>
            public Parameters()
            {
            }

            /// <summary>
            /// Constructor.
            /// </summary>
            /// <param name="f"></param>
            public Parameters(Form1 f)
            {
                Setup(f);
            }

            /// <summary>
            /// Fixed value heater percent.
            /// </summary>
            public double HeatPercent { get; set;}

            /// <summary>
            /// Watts.
            /// </summary>
            public double HeaterWatts { get; set; }

            /// <summary>
            /// Proportional
            /// </summary>
            public double P { get; set; }

            /// <summary>
            /// Integral
            /// </summary>
            public double I { get; set; }

            /// <summary>
            /// Derivative
            /// </summary>
            public double D { get; set; }

            /// <summary>
            /// Water heater time constant.
            /// </summary>
            public double TWH { get; set; }

            /// <summary>
            /// Water sensor time constant.
            /// </summary>
            public double TWS { get; set; }

            /// <summary>
            /// Water container time constant.
            /// </summary>
            public double TWC { get; set; }

            /// <summary>
            /// Container air time constant.
            /// </summary>
            public double TCA { get; set; }

            /// <summary>
            /// Water air time constant.
            /// </summary>
            public double TWA { get; set; }

            /// <summary>
            /// Water to food time constant.
            /// </summary>
            public double TWF { get; set; }

            /// <summary>
            /// Water start temperature.
            /// </summary>
            public double WaterStart { get; set; }

            /// <summary>
            /// Air start temperature.
            /// </summary>
            public double Air { get; set; }

            /// <summary>
            /// Water radius
            /// </summary>
            public double WaterRadius { get; set; }

            /// <summary>
            /// Water depth
            /// </summary>
            public double WaterDepth { get; set; }

            /// <summary>
            /// Duration
            /// </summary>
            public double Duration { get; set; }

            /// <summary>
            /// Hz
            /// </summary>
            public double Hz { get; set; }

            /// <summary>
            /// Control rate Hz.
            /// </summary>
            public double ControlRate { get; set; }

            /// <summary>
            /// Goal temperature
            /// </summary>
            public double Goal { get; set; }

            /// <summary>
            /// Control is binary or proportional
            /// </summary>
            public bool Binary { get; set; }

            /// <summary>
            /// Setup.
            /// </summary>
            /// <param name="f"></param>
            public void Setup(Form1 f)
            {
                HeatPercent = (double)f.numericUpDown11.Value;
                HeaterWatts = (double)f.numericUpDown12.Value;

                P = (double)f.numericUpDown9.Value;
                I = (double)f.numericUpDown10.Value;
                D = (double)f.numericUpDown13.Value;
                Goal = (double)f.numericUpDown16.Value;

                TWH = (double)f.numericUpDown7.Value;
                TWS = (double)f.numericUpDown8.Value;
                TWC = (double)f.numericUpDown6.Value;
                TCA = (double)f.numericUpDown5.Value;
                TWA = (double)f.numericUpDown17.Value;
                TWF = (double)f.numericUpDown19.Value;

                WaterStart = (double)f.numericUpDown1.Value;
                Air = (double)f.numericUpDown2.Value;
                WaterRadius = (double)f.numericUpDown3.Value;
                WaterDepth = (double)f.numericUpDown4.Value;
                Duration = (double)f.numericUpDown14.Value;
                Hz = (double)f.numericUpDown15.Value;
                ControlRate = (double)f.numericUpDown18.Value;
                Binary = !f.checkBox1.Checked;
            }
        }

        /// <summary>
        /// Constant heat.
        /// </summary>
        /// <param name="sender"></param>
        /// <param name="e"></param>
        private void button2_Click(object sender, EventArgs e)
        {
            RunSimulation(false);
        }


        /// <summary>
        /// PID test.
        /// </summary>
        /// <param name="sender"></param>
        /// <param name="e"></param>
        private void button1_Click(object sender, EventArgs e)
        {
            RunSimulation(true);
        }
    }
}

By viewing downloads associated with this article you agree to the Terms of Service and the article's licence.

If a file you wish to view isn't highlighted, and is a text file (not binary), please let us know and we'll add colourisation support for it.

License

This article, along with any associated source code and files, is licensed under The Code Project Open License (CPOL)


Written By
Technical Lead
United States United States
Phil is a Principal Software developer focusing on weird yet practical algorithms that run the gamut of embedded and desktop (PID loops, Kalman filters, FFTs, client-server SOAP bindings, ASIC design, communication protocols, game engines, robotics).

In his personal life he is a part time mad scientist, full time dad, and studies small circle jujitsu, plays guitar and piano.

Comments and Discussions