Differential equations are a cornerstone of modern mathematics. From quantum mechanics to population dynamics and stock market predictions, they play a crucial role in understanding the world around us. For this reason, courses on differential equations are core for many undergraduate degrees in the natural sciences, engineering and other fields.

Today, I am proud to announce a free interactive course, Introduction to Differential Equations, that will help students all over the world learn and master this difficult but tremendously useful branch of mathematics. The course will investigate the methods to solve a wide class of differential equations while making extensive use of the powerful functions to solve them in the Wolfram Language.

Want to get started? Explore the interactive course by clicking the following image before reading the rest of this blog post.

Motivation from History

Many of the great mathematicians in history have contributed to the development of the subject of differential equations. This began with Isaac Newton and Gottfried Wilhelm Leibniz and their introduction of the derivative in the seventeenth century. With the concept of the derivative of a function established, mathematicians could begin investigating equations involving functions and their derivatives. Such equations are known today as differential equations.

The study of differential equations brought about significant advancements in the eighteenth and nineteenth centuries. The work of Leonhard Euler formally established the solution methods for many types of differential equations during this time. Later, Pierre-Simon Laplace introduced his powerful integral transform, which can be applied to a wide range of differential equations. These methods of Euler, Laplace and others will be discussed throughout this course.

Overview

Students taking this course will receive a thorough introduction to differential equations. This includes standard topics like first- and second-order differential equations and initial value problems. Introduction to Differential Equations will also discuss some more advanced topics like special functions, the Laplace transform and systems of differential equations.

Here’s a sneak peek at some of the course topics (shown in the left-hand column):

I worked diligently to keep the course at a manageable length, and you will be able to finish watching all 31 videos and complete the eight short quizzes in about 10 hours.

Finally, I assume that students are familiar with college-level algebra, trigonometry and basic calculus. A background in linear algebra will be useful but not necessary, as any concepts from linear algebra that are needed will be reviewed in this course.

The next few sections of this blog post describe the different course components in greater detail.

Lessons

The backbone of the course is a set of 31 lessons, beginning with “What Is a Differential Equation?”. This introductory lesson—see an overview in the following video—includes a discussion about the applications of differential equations. This is followed by a brief history of differential equations and a discussion about where differential equations arise in the modern world, and concludes with an outline of the course.

Other lessons begin with a topic overview, such as direction fields or power series solutions, followed by a discussion of the key concepts with examples interspersed that illustrate the main ideas using Wolfram Language functions to verify solutions or visualize results.

The videos range from 10 to 15 minutes in length, and each video is accompanied by a transcript notebook displayed on the right-hand side of the screen. Copy and paste Wolfram Language input directly from the transcript notebook to the embedded scratch notebook to try the examples for yourself.

Exercises

Each lesson includes a set of five exercises to review the concepts covered in it. This course is designed for independent study, so detailed solutions are also provided for all exercises, as per this example:

The notebooks with the exercises are interactive, so students can try variations of each problem in the Wolfram Cloud. In particular, they are encouraged to change the coefficients of the differential equations, add terms and experience the awesome power of the Wolfram Language for solving differential equations.

Problem Sessions

Each section of the course includes a problem session that solves in detail several example problems that are related to the ideas developed in the section. These problem sessions are meant to reinforce the material that has been learned so students can solidify their understanding.

The following is an excerpt of the video for Problem Session 1:

Quizzes

Each course section concludes with a short, multiple-choice quiz with five problems. The quiz problems are at roughly the same level as those covered in the lessons, and a student who reviews the section thoroughly should have no difficulty in doing well on the quiz.

Students receive instant feedback about their quiz question responses, and they are encouraged to try any method—including hand or computer calculations—to solve them.

Course Certificate

Students should watch all the lessons and problem sessions and attempt the quizzes in the recommended sequence because course topics often rely on earlier concepts and techniques. At the end of the course, you can request a certificate of completion. A course certificate is earned after watching all the lessons and passing all the quizzes. It represents proficiency in the fundamentals of differential equations and adds value to your resume or social media profile.

An optional final exam is also available at the end of the course. Passing the exam allows you to receive an advanced level of certification.

A Building Block for Success

Mastering the fundamental concepts of differential equations is integral for students in the fields of science, technology, engineering and math (STEM) and is an immense intellectual achievement for anyone in any field. I hope that Introduction to Differential Equations will help you to achieve this mastery and contribute to your success in your chosen field. I enjoyed creating and teaching this course and welcome any comments regarding the current course as well as suggestions for future courses.

Acknowledgements

I would like to thank Devendra Kapadia, Cassidy Hinkle, Aaron Crenshaw, Amruta Behera, Joyce Tracewell, Veronica Mullen and Bob Owens for their dedicated work on various aspects (lessons, exercises and videos) of the course.

In a relatively popular marketing device in the past decade, chickens found their way into online advertising and TV commercials, where their impressive focusing and stabilization skills were displayed. Hold them gently and then move them up, down or rotate them slightly, and their eyesight stays at a constant level.
How do we build a device that mirrors those features? To find out, let’s embark on a journey from car suspensions to camera stabilization systems that reaches for those standards of excellence in engineering that can only be described as “majestic as the chicken.”
A Chicken Is Like a Spring
The first challenge we face is modeling such a system. An engineer could write a set of equations and plug them into a solver, but even for a relatively small model, this can produce a lot of results that are hard to share and visualize with peers because you can easily end with unconventional variable names, parameter choices and unit conventions.
Explore the contents of this article with a free Wolfram System Modeler trial.Given the choice, we would rather go to Wolfram System Modeler and use physical component–based modeling. If we want this model to have, say, a spring, then we drag a spring component and drop it on our canvas. We don t have to remember equations from college, and we get the added benefit of built-in animations.
It is not by accident that we mentioned springs. When we see a chicken’s neck expand and contract to stabilize its eyesight, it’s the first thing that comes to mind. A physicist would say that in a spring, energy that comes from motion is stored in the spring itself, but some could also be lost, depending on the amount of elasticity and dampening the system has.
This observation is what allows us to model structures that we parameterize in terms of how they handle incoming input energy, with part of that energy being stored and a part being lost. With this approach, we can build more complicated objects by putting together subsystems with different types of responses to inputs. Adding just the right elastic and dampening parameters can take us far on the path to modeling and controlling real systems.
Assembling the Model
Let’s make a simple model for the disturbances and dynamics involved in the stabilization of vertical motion. The system we use to mimic the response of the chicken is a camera mounted on top of a car. The disturbances caused by the road will transmit energy to the car and the camera. The objective for us is to design a control system that will negate those disturbances and allow us to have a smooth response.
We start by assembling a model of the system. We can emulate the interaction between the road, the car and the camera by connecting several spring-mass systems in series, each on top of the other. With the camera at the top, we let this system be affected by gravity, and we add two inputs. One represents the disturbance force, which acts from the bottom. The other one is a control input acting on the camera that will allow us to represent the device that stabilizes our recording setup.
(Download the model we created so you can evaluate the code in this post. Just make sure you save it to the same directory as the Wolfram Notebook you re using.)
Analyzing the Model
How do we know if we have made good choices for our parameter values? Let’s run some simulations and check that the results make sense. We start by importing the model:
#10005
If we ask for input variables, we get the disturbance force fw and a control input fa as inputs:
#10005
To improve camera performance, we need either a person or a device to hold the camera steady. We can skip the former because no human being is a match for a chicken’s reflexes. That’s the reason we have the actuator fa that is placed between the camera and the car body.
While controlling the camera position, we need some measurements that will help us determine how much force we need to apply to get the camera to its desired place. Two potentiometers can be a good fit for this job: one between the car tire and the body and the other one between the car body and the camera. These are the two outputs, ycb and ywb, of the model:
#10005
The state variables are the deviations of the car (xbe), wheels (xwe) and camera (xce) with respect to their equilibrium positions and their velocities:
#10005
We can see all these variables again when computing a linear state-space representation of the model:
#10005
If we give the system a sudden bump of 2000 N, our uncontrolled system oscillates very badly:
#10005
The camera jumped more than two centimeters!
Designing the Controller
By mimicking the mighty chicken, we got a system composed of a car and a camera. Its performance is way worse, however, than the chicken’s. It is time to control the camera position so that we can get better shots.
What control design can we implement? We can start by studying the poles of the system. For a linear system, these can be found by taking the state matrix and computing its eigenvalues, which will be complex numbers. The poles of the system play a crucial role because they have a big say in the characteristics of the system’s response.
It can be roughly claimed that the poles should always be on the left side of the complex plane. This condition guarantees that the structure will be stable. The further to the left the poles are, the better stability we have. We do not know precisely, however, which poles should be chosen. In other words, where should we draw a line, and which criteria should we apply?
Let’s consider the control effort. According to the chosen poles, the camera stabilizer force varies. If we want to have a quicker response, we need to aim for poles further left in the s plane. But this costs us more control effort. Thus, we have a tradeoff between the performance and the control effort.
The linear quadratic (LQ) controller technique helps us to express this tradeoff mathematically as an optimization problem that minimizes a cost function of the form . Here, x(t) and u(t) are state and input variables, and q and r are weighting matrices that capture the tradeoff. A small q matrix or large r matrix penalizes the states and will result in a controller that will produce a quick response at the expense of a larger control effort. On the other hand, a large q or r matrix penalizes the inputs and will result in a controller with a slower response and a smaller control effort.
The solution of the optimization problem is a simple linear equation of the form u *(t)= –k.x(t). Here, u *(t) is the optimal value of u(t) and k is a matrix of appropriate dimensions.
We have a built-in function, LQRegulatorGains , that has been recently updated to directly handle SystemModel objects and to automatically assemble the closed-loop system.
But how do we go about assigning values for the q and r matrices for this system? We can try and attenuate the whole system, including the vibrations in the tires and the car body. This would result, however, in an exorbitant control effort. Instead, let s focus on attenuating the vibrations of just the camera.
The camera s position (xce) and velocity states (vce) are the second and third states:
#10005
We will assign large values for the weights corresponding to these two states and have the values for the other states and the control input at a much lower value:
#10005
Now, compute a controller using these weights:
#10005
We asked for the data object, so LQRegulatorGains returns a SystemsModelControllerData object. This object can be used to query or compute several properties of the designed controller and closed-loop system:
#10005
Let’s look at the computed controller, which is typically of the form –k.x:
#10005
In accordance with the specified weights, we see that the second and third states, which are the camera s position and velocity, have the most influence on the control effort, while the third and sixth, which are the wheel s states, have the least effect.
The closed-loop model of the original model with the controller can be easily computed using ConnectSystemModelController :
#10005
We can then simulate it:
#10005
Take a look at its responses:
#10005
From the previous plot, it’s difficult to make out the camera’s response, which appears flat. So we’ll plot just that:
#10005
The camera on the uncontrolled system originally had 2 cm oscillations, and now that has been reduced by a factor of 104!
Let’s check the control effort as well:
#10005
It is less than 4 N, and it is a practical and cost-efficient control.
The weights q and r that were used previously were arrived at after a series of iterations. As we tweaked the weights, we evaluated the system s performance and the control effort. The notebook is quite a handy interface to carry out these explorations and record the values that look promising for further refinement and analysis.
Analyzing the Controller
Let’s analyze what the controller did to the system.
The uncontrolled system has a very oscillatory mode:
#10005
Now let’s look at the modes of the controlled system:
#10005
#10005
Notice that we have excluded the leftmost pole in the last plot to zoom in on the more interesting remaining ones. Two of the modes are essentially untouched, especially the seemingly problematic one that is closest to the imaginary axis. What the LQ controller did was to attenuate that mode’s response but did not remove the oscillations. That’s why we see those small oscillations up to 10 seconds, similar to the uncontrolled system.
What happens if we try to dampen that dominant mode a bit faster?
To do this, we are going to use the pole-placement approach to compute the previous gain matrix k. This is a single-input system, so the solution to the controller problem is unique and both the pole-placement and optimal LQ techniques will give the same result.
These are the poles of the closed-loop system that will not be changing:
#10005
We’ll create a function that will take the new eigenvalue, join it with the rest, compute the new controller and plot the new system responses and controller effort:
#10005
If we leave the dominant mode unchanged, we should get the same response as the LQ design:
#10005
If we try to move this mode ~0.02 rad/s, the control effort increases by a large factor:
#10005
If we move 0.2 rad/s, the control effort increases even more:
#10005
What we observe is that the oscillations in the controlled camera, although attenuated significantly, cannot be damped out without significantly increasing the control effort. The LQ design achieved a good balance, however, between attenuating the oscillations and keeping the control effort feasible.
Designing the Estimator
The controller we designed depends on the values of the positions of the three masses and their velocities. Using six measurements directly would not only require the same number of sensor readings but would also not compensate for any errors or noise in the readings. So we re going to use an estimator to provide the required values to the controller.
An estimator uses a set of measurements and the model of the system to estimate the values of the states. For example, we can estimate the number of people in an office building at any given time by counting the number of cars in its parking lot.
In our case, the two measurements are the relative positions of the camera and tire with respect to the car body:
#10005
ObservableModelQ tells us that it s possible to estimate the states from these two measurements and the model:
#10005
We are going to use the pole-placement approach to design the estimator. This works just as in the case of the controller, except that we’re placing the estimator’s poles. If we place the estimator’s poles well to the left of that of the controller’s, the errors between the actual and estimated values will decay fast. The following are the closed-loop poles:
#10005
After a few trials, we ended up placing the poles of the estimator at the following locations:
#10005
We then use the function EstimatorGains , which implements the pole-placement algorithm and computes the observer’s gains:
#10005
The estimator itself is a dynamical system, and the system is connected to the controller and estimator as shown here:
#10005
We don t have to do any assembling manually. Instead, we can use the function EstimatorRegulator to do that—and much more!—for us:
#10005
We can query the data object for the estimator model if needed:
#10005
But what we are really after is the complete closed-loop system:
#10005
Let’s now simulate the closed-loop system to see how well the system with the estimator and controller performs:
#10005
If we take a look at the system states, it appears to have responses similar to what we got without the estimator:
#10005
The control effort is also comparable to the controller without the estimator:
#10005
The response of the camera is now well damped, but it has a transient jump that quickly dies out:
#10005
Comparing it with the uncontrolled camera, we see that all the oscillations have been eliminated and the camera stabilized after a quick initial transient:
#10005
Improving Performance
Improving the performance of machine learning models often requires exploring the space of hyperparameters, metrics and data processing options. In a similar fashion, improving the performance of controllers requires testing and analysis. The key aspect when doing either is to have a framework where this investigation becomes natural and easy.
With an approach based on physical modeling in System Modeler, an expanding and continuously improving control systems functionality in the Wolfram Language, and their seamless integration, we have a framework of streamlined and powerful tools for modeling and controller design. And this is not the end of the line for them!
In this blog, we have shown these tools in action. We encourage you to try these out for the systems you want to investigate and design. For more information, explore the many domains covered in the Modelica Standard Library, our tech note on getting started with model simulation and analysis, our guide on system model analytics and design and our guide on control systems.
Check out Wolfram s System Modeler Modelica Library Store or sign up for a System Modeler free trial.
Engage with the code in this post by downloading the Wolfram Notebook
?Комментарии (0)

Explore the contents of this article with a free Wolfram System Modeler trial.You may not know what a Roland TR-808 Rhythm Composer is, but you have most certainly heard it. The TR-808 is a programmable drum machine released by Roland in 1980. The 808 is one of the most iconic drum machines and has been used in a wide variety of music, such as hip-hop, dance, soul, electro, pop and many more.
The 808 sounds like this:
#10005
Every 808 sound has its own personality. For over 40 years, it has fascinated musicians, who have used this drum machine in their music. If you are interested, you can find many videos about the 808, including this trailer for the 808 full-length documentary:
All the 808 sounds are generated using analog electronics. (Other drum machines may use samples, which are recordings of existing drums or percussion.) Over the years, I’ve been analyzing drum synthesis circuits, and the fact that all the drum sounds are generated by simple electronic components made the 808 very intriguing.
In this blog post, we are going to model and simulate the kick drum, one of the most emblematic drums of the 808, using Wolfram System Modeler. We’ll start by creating a component-level model to get an accurate reproduction of the behavior. (Download this zip folder for all the necessary files to experiment with the code in this post.) Later, we’ll simplify and optimize the model to compare the resulting sounds of the different versions against the 808 kick drum hardware.
Modeling the 808 Kick Drum
Here you can see the schematic I made for the bass drum:
#10005
(This schematic has small modifications compared to the original, but those modifications are not important for our analysis.)
When trying to understand a circuit like this, I start by looking for small subcircuits or patterns, such as filters, buffers, amplifiers or common operational amplifier circuits. When I’m not sure what the function of a circuit is, I simulate it and also corroborate the results against the real circuit. In this case, I have created my own replica of the circuit we’ll use as a reference:
#10005
Understanding Each Section of the Circuit
In the previous schematic, I have already divided the circuit into smaller subcircuits. We are going to analyze each of those sections by performing simulations individually, and at the end simulate the full circuit.
As we follow the signal flow, the first parts we will analyze are the processes of the Accent and Trigger signals. We’ll refer to this part of the circuit as Section 1.
#10005
Triggers are usually small pulses of around 1 ms that have an amplitude of 5 V. Accent signals define the velocity (intensity) of the drum sound hit. In this drum machine, there are three levels of intensity: a normal hit, an accent note and a ghost note. The voltage of the accent note defines the intensity. A large voltage (e.g. 14 V) will produce an accent note, while 8 V will be a normal note. In a similar fashion, a low voltage (e.g. 4 V) will produce a soft (ghost) note.
The Section 1 circuit consists of two transistors and a few resistors. We can simulate this circuit in System Modeler by using the generic transistor models included in the Modelica Standard Library because, in this specific case, the transistors act as switches. You can see in the following diagram the equivalent circuit made in System Modeler:
#10005
In this simulation, the Accent and Trigger signals are generated using the Ramp and Pulse blocks from the Modelica Library. The Accent goes from 0 V to 8 V, and the Trigger produces pulses of 5 V. We can see from the simulation results that when the Accent is above (approximately) 2 V, we get a pulse similar to the Trigger as output but with the amplitude of the Accent signal. This circuit behaves as a controlled switch where the Trigger signal controls the flow of the Accent signal:
#10005
We can see in the simulation that the switch is not perfectly linear. The accent needs to be at least ~2 V to get an output.
The next section of the circuit (Section 2) takes as input the pulse produced by the circuit we analyzed previously.
#10005
This circuit looks like a passive highpass filter, with the exception of the diode. The output of this stage is connected to the positive terminal of the op-amp. Thanks to that, we can assume it is connected to a high-impedance load.
We can recreate this circuit in System Modeler and use a generic diode model:
#10005
If we simulate this circuit using as input a 5 V pulse, we obtain the following results:
#10005
When the input pulse jumps from 0 V to 5 V, the output produces a peak that is almost as high as the input. When the input pulse returns to 0 V, the output has a small negative peak. By removing the diode and rerunning the simulation, we can see the purpose of the diode:
#10005
We can see a large negative peak in the output. The diode was put there to suppress the negative peak of this signal. The resulting pulse is now sent to the main resonator circuit (Section 3). This is the central part of the drum, and we’re going to analyze it next.
#10005
This part of the circuit is know as a bridged T-network or multiple feedback bandpass filter.
Let’s start by creating a simplified version of this circuit to analyze it:
#10005
In the previous circuit, we eliminated two of the inputs (the top of the original circuit), and the resistors of 47 k and 6.8 k have been combined into a single resistor of 53.8 k. We are using a small pulse of about 1 ms to excite (ping) the filter. The simulation results are displayed in the following plot:
#10005
We can see from the results that introducing a small pulse in the filter makes it oscillate (resonate) for a short time. The frequency of this oscillation is defined by the central frequency of the bandpass filter, which is defined by the values of the capacitors and resistors C1, C2 and R2 in our circuit.
You may remember that in this simulation, R2 in the original circuit is comprised of two resistors (47 k and 6.8 k), as shown in this circuit:
#10005
You may also notice there is a transistor connected in parallel with the resistor of 47 k (Section 6). This transistor seems to be operating as a switch that is used to short-circuit the resistor of 47 k. As mentioned before, changing this resistor value will make the filter resonate at a different frequency. In the following plot, we can see the results of short-circuiting the 47 k resistor, which changes the effective resistance from 53.8 k to just 6.8 k:
#10005
With a resistor of 6.8 k (47 k in short-circuit), the filter resonates at a higher frequency. A signal controls the switch (transistor in Section 6), but we are going to take a look at it later. Just remember that this signal is used to change the frequency of the filter.
At this stage, we can already listen to the drum sound by simulating our model (with the 6.8 k resistor). Note that you may need headphones or a decent sound system to listen to the low frequencies:
#10005
After the resonator, there is an operational amplifier that includes a potentiometer labeled as “Decay”. This is our Section 4, which is shown in this circuit:
#10005
In the drum context, “decay” refers to the amount of time the drum produces sound. A short decay sound is a fast sound, for example the sound produced by knocking a door. A long decay sound would be like strumming the string of a guitar, which usually takes a few seconds to stop making sound. Based on the potentiometer description, we can assume that this stage is in charge of controlling the decay of the drum.
This section is comprised of an operational amplifier working on inverting mode. If we ignore for a moment the 33 uF capacitor, we can see that we have the 500 k potentiometer connected in parallel with a 47 k resistor. The equivalent resistance of these two can be calculated with this formula:
#10005
In our case, the resistance of the potentiometer varies from 0 k to 500 k ohms depending on its position. Based on the previous formula, we can calculate the equivalent resistance in function of the position of the potentiometer (pos):
#10005
If we plot the equivalent resistance in function of the potentiometer position, we can see that this behaves as a logarithmic potentiometer, where its value goes from 0 k to ~42 k ohms:
#10005
As mentioned earlier, the operational amplifier is connected as an inverting amplifier. The formula for the gain of an inverting amplifier is given by the ratio of the feedback resistance over the input resistor. In our case, the feedback resistor is the equivalent resistance (Req) and the input resistor is 47 k; therefore, the DC gain of Section 4 is as follows:
#10005
Let’s perform a few simulations to corroborate the previous analysis and to also find out the purpose of the 33 uF capacitor in the feedback. In the following figure, you can see a slightly simplified version of Section 4. The only change is that we have dismissed the resistor connected to the positive input of the op-amp and fixed the position of the potentiometer to 10% (0.1):
#10005
Because this model has one input and one output, we can easily linearize it and create a Bode plot to analyze the frequency response:
#10005
From these results, we can observe that for very low frequencies (less than 0.1 rad/s), this circuit has a gain of 0 dB. For frequencies larger than 1 rad/s, the gain is given by the formula we calculated before. Basically, the 33 uF capacitor is there to let the DC pass unaffected, but attenuates all other frequencies according to the position of the potentiometer.
The output of Section 4 is connected (as feedback) to Section 3. This section is shown in the following figure. Let’s make a simulation model and test the effect of the decay:
#10005
In the next figure, you can see the simplified simulation model in which I have removed the 33 uF capacitor in Section 4. The decay potentiometer is set to 90%:
#10005
Now you can see the simulation results. In this case, the system resonates for a longer time and the next hit is triggered even before the drum stops producing sound:
#10005
By listening to the resulting sound, you can hear that it indeed has a longer decay:
#10005
So far, we have modeled a part of the main filter (resonator). To continue, we will delve into Section 5:
#10005
Section 5 uses two transistors that seem to be operating as switches. The first one takes the output of Section 1 (the Trigger signal). When the switch is closed, it quickly discharges the 100 nF capacitor. When the switch is open, the capacitor is charged through the 1 M resistor using the VCC voltage (around +15 V). The voltage of the capacitor controls the second transistor. When the voltage in the capacitor is low (discharged), the collector of the transistor (pin 1) is set to VCC through the 22 k resistor. As soon as the voltage of the capacitor polarizes the base of the second transistor, the lower terminal of the 22 k resistor is connected to ground.
Let’s simulate the circuit to have a better view of what is happening:
#10005
In this simulation, we are using as a trigger a pulse of 5 V with a duration of 1 ms. In the simulation results, we see how the capacitor discharges when the pulse is high. The capacitor takes around 6.3 ms before it closes the transistor T2, bringing the output to 0 V. In summary, Section 5 produces a pulse that is about 5.3 ms larger than the Trigger input. This signal is used by Sections 6 and 7:
#10005
The drum machine service manual mentions that this pulse is intended to be 4 ms, but in our simulation we get about 5.3 ms. This difference could be because we are taking a generic NPN transistor model that does not match the original transistor used. It is also important to consider that the real capacitors, transistor and resistors have variability, and when combining these three components, the timing is most likely to be off. In our simulation model, we can easily match the timing by reducing the 1 M ohm resistor.
We have seen previously that the transistor in Section 6 bypasses the 47 k resistor of the main resonator, making it oscillate faster. Because the transistor is controlled by the pulse produced by Section 5, we can conclude that the drum oscillates faster for a short period and then returns to its main oscillation frequency:
#10005
That change of pitch in the resonator is there to emulate what happens in a real kick drum. When hitting a real kick drum with a mallet, the membrane of the drum is stretched, which increases the tension, making the pitch higher for a very short period of time. After a few vibrations, the drum returns to the original tension and continues vibrating until the energy is dissipated.
Section 7 receives the pulse of Section 5 too. Let’s create a simplified circuit and simulate it to understand what’s happening:
#10005
In the previous circuit, we recreated Section 7 plus an emulation of the pulse produced by Section 5. In the following plot, you can see the simulation results. This circuit produces a negative pulse with a similar width as the input pulse. The way this circuit works is that when the transistor T1 is open, the capacitor C1 is charged through R1 and D1 to +15 V, and when the transistor is closed, it connects the positive pin of the capacitor to ground. This applies the inverted voltage of the capacitor (–15 V) into D1 and the output. This signal is used to excite the resonator in addition to the trigger pulse:
#10005
The last stage we are going to analyze is Section 8. We are not going to cover Section 9 because it is very simple and is basically a level control and DC offset removal. Section 8 has a potentiometer labeled “Tone”. It is easy to see that this is simply a lowpass RC filter. One interesting thing to notice is that the potentiometer is connected in parallel with a 10 k resistor. This has a similar behavior as the Decay control we analyzed before. It makes it a logarithmic-like potentiometer, but the curve is a little bit different. In this case, the equivalent resistance is given by the following equation:
#10005
The relation between the equivalent resistance and the position of the knob looks as follows. Using that value, we can easily calculate the cutoff frequency of the RC filter, which is given by the formula:
#10005
If we substitute the values of the components, we get the following value of the cutoff frequency in function of the potentiometer position:
#10005
Now that we have a better understanding of the different stages of the circuit, let’s put them all together to simulate the full circuit:
#10005
In this simulation, we have set the following parameters: Tone 10%, Decay 35% and Accent 8 V. As mentioned during the analysis of Section 5, we have reduced the 1 M resistor to 900 k to get a pulse closer to 4 ms. We can see the result in the following plot:
#10005
We can see that after the trigger, the resonator is excited by a short pulse (Section 2) whose amplitude is given by the Accent amplitude. This makes the resonator start oscillating. The 4 ms pulse of Section 5 bypasses the 47 k resistor of the resonator by short-circuiting the transistor of Section 6. The capacitor of Section 7 starts charging to +15 V. When the resonator is about to complete half a cycle (which is around 4 ms), the pulse of Section 5 goes low. This introduces in the resonator a pulse of –15 V (Section 7 output), and at the same time, the frequency of the resonator changes when reintroducing the 47 k resistor.
In summary, we can see from the output that there is a small positive “click” followed by a negative half-cycle of a sine wave. This half-cycle is about twice the base frequency of the resonator. After that half-cycle, the resonator oscillates at the base frequency until the energy dissipates according to the set Decay.
We can listen to the sound produced by this simulation in the following audio output. (Again, because it is a low frequency, you may need to play it using a decent sound system.)
#10005
Here are a few insights about this simulation model:
It contains 160 equations.
It requires solving 7 linear systems of equations and 1 nonlinear system.
It simulates 7 semiconductor components (5 transistors and 2 diodes).
On my computer, it takes 3.5 s to run a simulation of 4 s.
In our previous analysis, we discovered that practically all semiconductors in the circuit are used as switches. We could create a more efficient model by simplifying all the circuit sections according to their behavior. But before doing that, let us take a look at how this model matches our real kick drum.
Comparing to the Real Circuit
In the beginning of this post, we saw a picture of the circuit I built to use as a reference. There are a few things to consider when comparing the model against the real circuit.
First, my circuit uses resistors with a tolerance of 1%, while some of the capacitors have 5% and others 1%. My replica of the drum is powered by ±12 V coming from my Eurorack power supply, while the original drum machine used ±15 V. This difference may have a small impact on Sections 5 and 7, according to our analysis. But these differences can be compensated either in the simulation model or in the real circuit. Finally, I used general-purpose transistors that do not match the ones used in the original drum machine, but because they are used as switches, I expect the differences to be minimal. It is important to note that in the real circuit, it is hard to set the exact values of Tone and Decay, and this may have a small impact on the measurements.
In the following image, you can see the measurements taken from the input Trigger (blue) and the pulse of Section 5 (yellow). You can see that the total length of the pulse is around 6.2 ms. This number is very close to the 6.3 ms of our simulation. These numbers are close, however, to each other for different reasons. As mentioned before, this circuit is powered with ±12 V, which will make the capacitor of Section 5 charge slower, therefore making the pulse longer. The transistor model needs a slightly larger voltage to turn on; therefore, the pulse is longer:
#10005
In the next screenshot, we can see the output of the drum (blue). We can see in this measurement that during the 6.2 ms, the resonator almost completes half a cycle. The 6.2 ms end when the resonator is close to the positive peak. This implies that the resonator is oscillating faster than the simulated resonator:
#10005
In the following plot, we can see the measured drum sound (blue) and the simulated sound (orange) side by side:
#10005
In this plot, it is easier to see that the real resonator oscillates faster than the simulated one. Now let’s listen to the sounds:
If you have a good ear, you may notice that, indeed, the real drum sounds slightly higher in pitch. Apart from that, they sound very similar.
At this point, there are a few paths that we can follow:
Tweaking the model so it matches the real drum.
Fixing the real drum so it matches the service manual (making the pulse 4 ms).
Making the model match the service manual by adjusting the frequency and fixing the 4 ms pulse.
For this example, I will try to make the model match the real drum by tweaking some of the parameters. After changing the Tone and Decay parameters and modifying the resistor values of the resonator, I got the following sound:
This tweaked model is a better match. There are still a few improvements that we could make to our model if we would like to get a better matching of the waveforms. For example, in the real circuit, the resonator frequency seems to change more gradually, so it is faster in the beginning until it settles in the main frequency. Our simulation model settles sooner. This could be caused by the nonlinearities of the transistors or the op-amp. This behavior could be simulated by tweaking the transition curve when switching the resistor in the resonator on or off. Based on the previous results, however, the sound of our model is close enough.
In the next section, we are going to try to simplify our model even further and also try to improve it based on the observed behavior.
Simplifying the Simulation Model
As already mentioned, the simulation model presented can be optimized by using simplified models of the components. We can also take advantage of the power of the Modelica language in System Modeler to create simplified behavioral models that match the behavior of the original circuit. The objective of this is to create a simpler and faster model that we can port to an audio platform and use to create actual music. In my previous blog post in this series, I showed how I create plugins for the VCV Rack open-source Eurorack simulator.
Because we already have the circuit split into sections, we can try to create an equivalent simplified model for each section. We’ll use causal connectors to communicate each of the sections of the circuit. This means instead of using electric terminals (acausal connectors), we will use inputs and outputs.
Let’s start with Section 1. When analyzing this section, we concluded that it behaves as a controlled switch. The Trigger signal, which can be represented as a Boolean signal, controls the flow of the Accent signal. Using System Modeler, we could replicate this behavior by using some of the existing components. In this case, however, I have decided to create my own model to show you how to make custom models using the Modelica language.
I start by creating a new model in System Modeler. (Check out the System Modeler tutorials to learn more about model creation.) Then I drag the inputs and outputs into the Icon View to locate them in my preferred position:
#10005
Once I have the icon for my model, I can switch to Text View to type the code of my behavioral model:
#10005
You can see in the code that this model consists of a single equation that defines the value of the output. You may notice that as soon as the Trigger signal becomes “true,” the Accent is sent to the output. This model differs from the original when the Accent signal is less than 2 V. In that case, the original model did not produce an output. This is not an issue because the Accent signal is supposed to have a range from 4 V to 14 V. On the contrary, this simple model would allow us to produce sounds with an Accent less than 2 V. In the following plot, you can see the results. Modeling it this way removes two transistors:
#10005
Continuing with Section 2, this model is easier to create by using existing components. If you recall, this section used a diode to suppress the negative peak produced by the highpass filter. We can create an approximated behavior by removing the diode and adding a limiter at the output:
#10005
Compared to the original model, this simplified version completely blocks the negative output. The biggest difference is that the capacitor is not being discharged by the diode. Therefore, if a new trigger comes before the capacitor is discharged, the output will be slightly different. This is not a big issue since we do not expect to have two triggers in less than 300 us, which is the time the capacitor takes to discharge:
#10005
By removing this diode, we have removed from the complete model the only system that required a nonlinear solver when simulating.
Before going into Section 3, we are going to continue to Section 5.
Section 5 is in charge of producing a 4–6 ms pulse for every Trigger. We can easily describe this behavior with a snippet of Modelica code:
#10005
This simple model produces the pulse and outputs it in two ways: one as a Boolean signal and the other as a Real. In this specific case, the pulse is 5.3 ms to match our previous simulation model. By using this simplified model, we get rid of two more transistors.
Section 7 is modeled using a combined graphical and textual model. In this graphical model, you can see one capacitor and a single resistor. The resistor acts as a load and discharges the voltage of the capacitor:
#10005
The capacitor is charged to –11 V when the signal of Section 5 (the 4 ms pulse) ends. To describe that behavior, we use a “when” equation in Modelica. This equation reinitializes the state of the capacitor to the target voltage. Here is the code part of this model:
#10005
By making this change, we get rid of the second diode. In the next plot, you can see the results of the behavioral models we created for Sections 5 and 7:
#10005
Now let’s go back to Section 3. If you take a look at the following diagram, we have combined Sections 3, 4 and 6:
#10005
The transistor used in Section 6 to change the oscillation frequency has been replaced by a variable resistor that takes the Section 6 signal and uses it to define the corresponding resistor value. The capacitor in Section 6 that was used to let the DC pass unaffected is removed.
Finally, Section 8 is very similar to what we had before: a simple lowpass RC filter. You can see it in this diagram:
#10005
When we put all the simplified models together, we end up with the following simplified model:
#10005
Before seeing the result of the new model, let’s take a look at some insights:
It contains 55 equations (compared to 160 in the original model).
It requires solving 3 linear systems and 0 nonlinear systems.
When simulating 4 s, it takes 0.35 s (almost 10x faster than the original).
This plot shows the results:
#10005
When looking at the output signal, it is hard to spot the differences. Listen to both sounds:
Both models sound extremely similar. In a blind test, it would be hard to determine which one is more accurate and which one is simplified. If we look at the simulation results, however, we can spot small differences. Here, the original model is in blue and the simplified in orange:
#10005
Here is a close-up of the initial transient:
#10005
I tweaked a few parameters in the simplified model to better match them, including the value of the main resistor of the resonator, and added a few gain blocks. Some of the differences in the results may come from the fact that we replaced all the transistors with ideal switching devices. We removed one capacitor in Section 4. In addition, all our sections are completely decoupled from each other, while in the original circuit, because of the nature of electric components, one section can affect the behavior of others though the current flow.
When we compared the original model against the real drum, we mentioned that the real resonator did not make an immediate transition to the base frequency. In the simplified model, we can replicate that behavior by changing the Section 3 resistor value gradually, for example by using a lowpass filter. In addition, I inserted a gentle nonlinearity at the output to limit the amplitude of the initial peak.
Listen to the real drum and improved model:
#10005
#10005
The following plot shows a comparison of the waveform produced by the model and the real drum:
#10005
You can see in the previous comparison that, overall, this model is a good match to the real drum. It still fails to reproduce, however, all the effects in the initial transient:
#10005
At the end, the selection of the model variations presented in this post depends completely on the application. If we want to run these models to consume as little CPU as possible, we can further simplify them. If we want to achieve full realism of the sound, it is at the expense of computational power.
In summary, the simplified model is a good approximation to the full circuit, and it can be tweaked to improve the matching with either the real hardware or the complete simulation model.
Extending the Basic Functionality
Based on our analysis, we were able to find a few components whose values may have a significant effect on the sound. We could replace these fixed components with variable versions to expand the sonic capabilities of the drum. For example, we could change the resistor setting the pitch of the drum by a potentiometer (see the next figure). That would allow us to change the pitch of the drum. A value of 50 k will give us approximately one octave of range:
#10005
Other possible modifications would be to change the length of the Section 5 pulse, but we are not going to cover it here.
The following is a demo of the sounds that could be produced if we changed all the parameters of the drum (including the pitch):
#10005
#10005
Optimizing the Model for Real-Time Execution
Implementing the simplified model is still a nontrivial task. There are some sections, like the main resonator, that we have described using electric components. As in my previous blog post, the model would be easier to implement in a language like C++ or Vult if it is in ODE form. We are going to see how this can be done for the main resonator. The other sections are easier, and a similar approach can be followed.
The component for the main resonator is made in such a way that it uses only linear components. This will make it easier to reduce the equations to a minimal form by using the functions provided in the Wolfram Language to manipulate System Modeler models:
#10005
#10005
#10005
#10005
These three equations are easily translated to the implementation language. In this blog post, we are not going to cover the implementation process, but you can refer to my previous post, “Digital Vintage Sound: Modeling Analog Synthesizers with the Wolfram Language and System Modeler.”
Additional Resources
The way the circuits behind drum modules work is fascinating. Over the last year, I have analyzed the internals of many of them. This work has resulted in a collaboration between myself and the creator of VCV Rack to release an audio plugin that simulates a full drum machine.
By using the modeling capabilities of System Modeler, I was able to create accurate models out of real drum modules that I used to understand the behavior of the circuits. Once I had a good understanding of the circuit, it was possible to create and test simplified models that behave close to the original circuit. The simplified models take advantage of the mixed graphical and textual modeling approach to easily and efficiently reproduce the original models.
You can even listen to some of the demos made with these drum modules in this YouTube video:
Learn more about modeling electronic circuits using System Modeler with Wolfram’s free CollegeDigitalElectronics download at the Library Store or sign up for a free System Modeler trial.
Engage with the code in this post by downloading the Wolfram Notebook
?Комментарии (0)

A version of this post was originally published on Tech-Based Teaching, a blog about computational thinking, educational technology and the spaces in between. Rather than prioritizing a single discipline, Tech-Based Teaching aims to show how edtech can cultivate learning for all students.
A blank page, bane to any writer!
If you’re trying to write a story, be it for National Novel Writing Month or just for fun, you’ll have to face a blank page eventually. The seeds of an idea can help your story grow, blooming into a sweet rose of romance or a carnivorous tale of horror. Without those ideas, all that’s left is a blinking cursor… and frustration.
If you’re stuck on what to write, maybe a little randomness can help. Using some of the Wolfram Language’s Random functions, you can brainstorm ideas, pull from image prompts, create mood boards and even chart the rise and fall of a plot.
Ready? Let’s write!
Getting Our Setup: Genre and Subject
What genre are we writing in? Genre is a way of framing ideas, featuring useful tropes we can use or deconstruct. To get a genre, let’s fill up a word list with some genres we enjoy writing. We can then pick one at random with RandomChoice .
Here is what we’d type for a choice among sci-fi, fantasy and horror, as well as the random result:
#10005
What is our story about? Let’s use ResourceFunction [" RandomWikipediaData " ][] to pull up a random Wikipedia page:
#10005
Hm. That’s okay, but if we reevaluate the function, maybe we’ll find something better:
#10005
Yes! Who doesn’t love a fjord?
Note that once we have the function, we can highlight the cell and evaluate it over and over until we get an inspiring result. (Or you could write on hard mode and use whatever result you get, even if it s the bio of an obscure Polish shoemaker.)
If you re not inspired through the RandomWikipediaData function, you can get some more ideas through a random list of words with ResourceFunction["RandomText"][3]:
#10005
Two friends travel to Greenland, and due to mixed signals, they have a falling out, resulting in sulkiness once their boat’s anchor malfunctions….
A sulky scientist sets up anchorage at a fjord and finds an unusual signal….
Once you’ve gotten an idea or two, you can find your story’s theme by changing the number at the end of the function to [1]. This will show a single, thematic word.
Let’s try it. Our theme is:
#10005
Not your typical theme, but perhaps it’s a metaphor for hidden depths. The sulky scientist isn’t so one-dimensional after all!
Visualizing the Story with Random Images
Having words in mind is good, but another way to keep the ideas flowing is to use random images. There are several random image functions available, so let’s try one and see what happens.
If we type in ResourceFunction["GetLoremFlickrImage"][], we can get a Creative Commons–licensed image:
#10005
If there’s an internet, there’s a cat. Guess the scientist has a feline friend. At least the fjord has plenty of fish.
Images can help with all sorts of visualizations. For example, if there are people in the images, they can act as friends and enemies (to lovers). If there’s a location, it can be a scene in our story. A random object can act as a MacGuffin or Chekhov’s gun to ramp up the story’s tension.
How about looking at the story as a whole? Using ImageAssemble [ Table [ResourceFunction[" RandomPhoto "][150,150],3,3]] gives us nine random photos arranged in a table to create a mood board. Mood boards, like theme and genre, reflect the story’s core ideas and reader expectations:
#10005
It’s a bit of work to connect the images, but there’s bleakness, prickliness, open water and looking up from the depths. There’s no horror just yet, but the fog is ominous.
Note that by priming the story with the previous random prompts, it’s easier to see where this mood board might connect with those initial seeds. This same set of random images could just as easily skew toward a rom-com where a surly, high-achieving businessman moves back to his hometown by the sea, only to discover a special someone who protects wild moose….
What a Character!
Now that we have an idea of where and what our story is, let’s figure out who the story is about.
Let’s name our characters with ResourceFunction[" RandomPetName "][], which will give us a random name. If we run the function in two different cells, we can get two random names at the same time:
#10005
#10005
Salem Adams, sulky scientist, and Flora, feline first officer, star in a tale of horror as they investigate a mysterious signal from the depths of Greenland’s fjords….
Random functions aren’t just good for naming your characters: they can describe what your characters do. Say you’re writing a romance. What’s your love interest’s job? Use ResourceFunction[" RandomPretentiousJobTitle "][] to find out:
#10005
C’mon, we can get more pretentious:
#10005
For all the sci-fi writers, let’s throw in an alien name for fun. Use ResourceFunction[" RandomString "][] to create a random string of letters. Put a number in the brackets for a specific number of letters:
#10005
His friends call him Ulb. His enemies never even learn his name….
Hey, Plot, Where Are You Going?
At this point, your fingers might be flying over the keys, drawn into this random story. But what’s the shape of your story?
There are many ways to look at narrative and story structure. There are five-act stories that have a rise and fall. There’s the monomyth, which cycles through patterns and archetypes. There are narratives that loop and reshape themselves around characterization, focusing on snatches of time.
We’re going to let Graphics [ BezierCurve [ RandomReal [1,{5,2}]]]—a function that draws a random B?zier curve through a certain number of points in a certain area of space—tell us the path of our plot. Maybe it’s the level of excitement in our narrative. Maybe it’s character growth over time. Either way, we can use the random path as inspiration.
Let’s see some examples:
#10005
The situation steadily worsens before history repeats itself, leaving Salem once more facing peril.
#10005
Salem is about to lose his grant funding, and he’s hoping that the fjord trip will result in something big. Too bad he gets his wish, ending up at a lower point than when he started.
#10005
It’s one step forward, one step back, but Salem is never happy.
Feel free to change the numbers and see what happens. If you want to dive deeper into random scribbles, check out this explainer that explores all the parts of the function and takes you into the 3D space. (How might that change your story’s plot structure? Perhaps converging narratives?)Комментарии (0)

A couple of years ago, I wrote a piece outlining why I think that open source isn’t the right business model for Wolfram’s core tech. It generated some (mostly reasonable) debate about the benefits of different models.

Reflecting on that debate recently, I realized that most of the expected practical benefits of open-source software are also strongly apparent in Wolfram tech despite our non–open source approach. So if you can forgive the clickbait title, I thought I’d walk through them:

Because our central business model is making great software and licensing it for money, it surprises many people to learn that there are quite a few ways to use the Wolfram Language for free. Of course, Wolfram|Alpha is free, but I mean the full Wolfram Language. The easiest way is to create a free Wolfram Cloud account. You have access to the whole language, through a browser or via APIs, for free. True, there are some memory, CPU time and storage limits as it costs us money when you use the free cloud, and you can pay to upgrade those.

Prefer local? Try the Raspberry Pi version. The complete Wolfram Language runs on the $5 computer for free with some commercial use restrictions. Want to run it on a PC? Then there is the Wolfram Engine for macOS, Windows or Linux, free during the development phase of your project, or Wolfram Player, free for running code only but not for writing new code.

So it’s not the “free to do anything” of open-source software, but it is free in many cases.

On top of that, one has to remember that the full tech stack is free at the point of use to millions of people thanks to their institutions’ support. Most good universities have site licenses, so any student or faculty can use the technology unrestricted without personal fees. We even have country arrangements—if you are a student, teacher or academic researcher in Egypt, you have free Mathematica—that’s about 40 million people right there. You might say, “That’s not free because someone is paying,” but isn’t that true for open source? Millions of dollars have been spent so far on Jupyter and related projects just to make thin versions of Mathematica notebooks, funded in part—annoyingly!—by my own personal taxes.

2. You Can View the Source Code

A large and increasing fraction of the Wolfram tech stack is written in the Wolfram Language, and all of that source code is viewable. Just turn off the ReadProtected attribute and ask for the definition. That’s been true since Version 1, but for some time, there was an undocumented but widely known internal function to allow you to browse the definitions in a point-and-click way. That function is now documented in the Wolfram Function Repository as ResourceFunction["PrintDefinitions"].

Here, for example, is the result of ResourceFunction["PrintDefinitions"][URLRead]:

✕

Click functions in the code to see their definitions. Older or core functions written in C are not viewable this way, but of course, some of those are calling actual open-source libraries like MXNet, MKL, GMP and ImageJ, all of whose sources are available to view elsewhere.

3. You Can Modify It

Any public or internal function in the Wolfram Language can be replaced. Look up the source, edit it and evaluate it to replace the built-in rule. More than that, however, the fundamental design of the Wolfram Language, with its universal support for operator overloading, allows you to replace or augment any functionality without even having to look at the source.

For example, suppose I think that I have a better implementation of Sin for small reals, but I don’t want to replace the built-in behavior for large reals, complexes or symbolic arguments:

✕

✕

✕

You can see that my new implementation is being used without having disturbed the other behavior. That’s a lot easier than editing source code:

✕

4. You Can Contribute

As I explained in my open-source post, we don’t really believe in user contributions to the core language—although we do like well-described bug reports. But there are many other ways to contribute to the Wolfram Language ecosystem. At the time of writing this post, there are more than two thousand additional Wolfram Language functions in the aforementioned Wolfram Function Repository, contributed by users as well as Wolfram developers. All of them are immediately available within the language as if they were built in (but in reality auto-downloaded and installed on first use) and all with source code.

For example, here is one that I contributed for creating geo surface plots with embedded maps:

✕

There are also more than 12 thousand Demonstrations that have been contributed for educational purposes, again with source code.

You don’t have to use our delivery systems. You can find many examples of open-source Wolfram Language packages shared in other ways, such as GitHub. Some can be seen at Wolfram Community and PackageData, and we are preparing a new way to share larger packages of Wolfram Language code.

5. There Is a Supportive (and Supported) Community

There are actually two organized communities, each with a slightly different emphasis.

The completely independent Mathematica & Wolfram Language site is at Stack Exchange. Structured in a Q&A format, it has handled more than 80 thousand questions. It’s the only Stack Exchange site dedicated to a single programming language and almost the only one dedicated to a single technology stack.

We also host and actively support Wolfram Community. With over 30 thousand members, it’s a more general site where people can share their work, discuss approaches and ask for help. It’s a great place to find collaborators for packages that you might want to contribute to the community.

You will find many contributions on both sites from Wolfram developers.

6. There Is Active Development

Aside from the minor patches that we regularly push out, we shipped six full releases in 2020 and 2021. The top four averaged about 170 completely new functions on top of many improvements to existing functionality. That number is more than the total operator count in the core Python language.

We Practice “Open Design”

This one wasn’t on the list of things that people said were good about open source because even though you might expect that to be the place where open design happens, it really isn’t. You can often see the history of bugs and pull requests, but not much discussion about design, strategy or why decisions are made.

We also broadcast design review meetings live. You can hear how we decide what functions get into the Wolfram Language, how we work out the minimal set with maximal functionality and how we make sure they fit together as a whole (all those things I talked about in my previous post). You can even ask questions live in the meetings or if you have a good idea for something we are stuck on, suggest it. Sometimes the meetings are full of elegant insights and sometimes they get quite heated, but often they just reveal how hard good design is, as we grind through some difficult, obscure—but important—detail. And you can watch the archive of over five hundred hours of livestreams.

So, while we are not about “free and open-source software” in the usual sense, Wolfram tech is sometimes free, mostly open and always well designed.

The Wolfram Language is renowned for simplicity and brevity, and nowhere was that more apparent than at the 10th annual One-Liner Competition, held at the Wolfram Virtual Technology Conference. The contest challenges conference attendees to create the best program possible in 128 characters or fewer (the original length limit of a tweet). With prizes awarded for the three best submissions, competition this year was fierce, but the judges, with only minor bloodshed, were able to settle on a slate of awardees.
As a reminder, submissions were to be linear strings of 128 characters or fewer, with no 2D typesetting constructs allowed. The judging criteria were highly personal and dependent on the capricious whims of the six individual judges (with impartial moderation provided by emcee Daniel Lichtblau), but it included brevity, aesthetics, originality and surprising use of the Wolfram Language. Wolfram employees were not eligible to compete.
Without further ado, here are the winners .
Honorable Mention
Benny Cheng: Bifurcation Diagram (126 characters)
This submission holds a special place in the heart of this judge, who fondly remembers producing the same diagram in teenage coding experiments in the 1990s using QBasic. What is not so fondly remembered is the length of time and code needed to make such drawings given the tools available to an amateur of the time. How far we’ve come!
#10005
The judging panel also agreed that the image output just looks cool. The use of 16 characters to add Background Black seems worth it. October is Spooky Month, after all.
Third Place
James Lu: Minute Repeater Computa UltraTerse (122 characters)
Our third-place winner takes the telling of a story to a higher level in the description of his submission. Here are his own words:
.product-details { background: #efefef; padding: 10px 20px 0;}
“Minute Repeater” is a combination of mechanical watches and one of the most fascinating accomplishments of Haute Horlogerie . At the press of a button, a minute repeater watch chimes out the hour, followed by the quarter hour and minutes. For example, if the time is 3:34am, then the minute repeater will sound three low tones representing three hours, two sequences of dual-tones representing 2 ? 15= 30 minutes, and four high tones representing four minutes: “ding, ding, ding; ding+dong, ding+dong; ding, ding, ding, ding.” Such watches originated before night illumination became widely available and were invaluable for the visually impaired. In the current day, while minute repeater watches have lost their practical utility, they are still highly sought after by watch enthusiasts due to their intricacy and novelty.
The mechanisms underlying minute repeater watches are ingenious and require both technical facility as well as musicality on the part of the watchmaker to listen to and confirm the harmony of the striking sounds. Furthermore, the pinnacle of achievement for minute repeater movements is their ability to fit into a small watch case. As a tribute to this mechanism, I’ve utilized the computational “gears” of the Wolfram Language and implemented a terse code succinct enough to fit into a 128-character “case.” At the press of “evaluate,” this minute repeater watch tells time in a manner that harkens back to an earlier age but relies on modern computation.
#10005
Judges loved the history lesson and appreciated the accessibility applications of this program for the sight-impaired. Note the clever use of infix notation in # 1~ Table ~#2 to save one character over Table[#1,#2].
Second Place
Michael Sollami: Poor Man s Vacation (126 characters)
Everyone needs a vacation now and then, right? And in the time of COVID, perhaps virtual travel is the best we can hope for. This submission allows the user to dynamically choose a vacation spot somewhere on the globe and then displays images obtained through WebImageSearch to show you what you might see if you traveled there. Note the definition d= Dynamic here, which allows an overall savings of four characters. By joining code with \[NewLine] characters rather than semicolons, as in a CompoundExpression , the author allows us to view multiple linked Dynamic expressions simultaneously—a neat idea!
#10005
First Place
Michael Sollami: Chessboard in a Tweet (128 characters)
Who would have thought that one could build a functional chessboard in the Wolfram Language with only 128 characters?! In this elegant code, first-place winner Michael Sollami does just that. Take a look at the code and its output:
#10005
Remarkably, one can play a game on this board. The click-and-drag functionality needed to interact with the pieces in the front end is provided entirely through Locator used within Graphics , showing the power of the Wolfram Language when coupled with the dynamic capabilities of the Wolfram front end. Judges were very impressed with the functional programming here to make everything fit within the character limit. The use of Raster and Array with an appropriate grayscale function condenses the code for visualizing the board itself (without pieces) to a mere 27 characters. In the first argument of Array, the author takes clever advantage of the prefix form of Plus , whereby +## is equivalent to #1+#2, but uses three fewer characters. Also of note is that the chess piece figures, being standard UTF-8 characters, count as one character apiece as per the contest rules. Why not grab a friend and try a game?
What Are You Waiting For?
After seeing what the Wolfram Language can do in a very small amount of code check out all of this year s submissions we hope you’ve been inspired to have a little fun and create something unique of your own. If you come up with your own one-liner, Wolfram Tweet-a-Program, described here, is a great way to share it with the world.
Thanks to all who participated, and congratulations to this year’s winners!
Engage with the code in this post by downloading the Wolfram Notebook
?Комментарии (2)

English bell ringing (called change ringing) has many connections to mathematics, notably to group theory and Hamiltonian cycles. My wife, Joan Hutchinson, is an ardent bell ringer (having rung in both England and North America), and I knew the basics of this ancient craft. A recent puzzle book by Mark Davies [1] inspired me to bring Mathematica’s integer-linear programming (ILP) capabilities to bear, but I wanted to go beyond puzzles and develop a new ringing method that would be of interest to the bell-ringing community.

My idea of finding a method with a large number of “wraps” was a success, and the method, which I called Wolf Wrap, has now been rung by a British band of ringers; this means it was formally entered on the Ringing World’s list of recognized methods. Here I will describe how one can interpret the problem in terms of a linear system of inequalities and so apply ILP, which is accessed via Mathematica’s LinearOptimization function.

1. Bell-Ringing Basics

The English hang their bells in a unique way: each bell is attached to a wheel that rotates when a rope is pulled by the ringer. The wheel bell in the tower swings through a full revolution from mouth upward to mouth upwardfor each sounded tone. There is one ringer for each bell in the tower; when methods are rung on hand bells, each ringer controls two bells.

This is the “ringing chamber” at Stoke Gabriel Parish Church, Devon, England. The ringer pulls on the rope to turn the wheel to which the bell is attached.

These are the bells of St Bees Priory shown in the up position. When rung, they swing through a full circle from mouth upward around to mouth upward, and then back again.

A typical tower has eight bells, ranging in pitch from the treble (bell #1) to the tenor (#8). When rung in the natural order 12345678 (called rounds), the tones form a descending major scale. Ringing methods require ringing the eight bells in the order corresponding to a permutation of {1, 2, 3, …, 8} and continuing through dozens to thousands of such “bars” of music (I will call each bar a row). The essential rules are:

The first row 12345678 is usually played twice to establish rhythm and parity, but will be given only once in the work here.

A row must never repeat a previous row, except that the final row is the same as the first. This property is called truth.

A bell can never move more than one position from its current position.

Rule 3 is a consequence of the physical nature of the wheels. It takes time for the wheel to rotate, and proper rhythm means that a bell cannot move far from one row to the next: it either stays in place (called making places) or moves one place to the right or left.

To see and hear more about how this works in a bell tower, watch these YouTube videos:

(This article has more information on the history of change ringing.)

Methods have been developed so that—in the case of seven bells, for example—a sequence of 5,040 permutations can be rung (from memory: no written notes), starting and ending with rounds and never otherwise repeating a row. Really there are 5,041 permutations, but the number of changes—the moves from one row to the next—is 5,040. This is called a full extent on seven bells, as all 7! permutations appear. A full extent on eight bells has 40,320 changes, too many to be practical.

Any sequence of changes numbering at least 5,000 is called a peal (a quarter-peal has at least 1,250 changes). Ringing a peal or quarter-peal is a significant mental and physical challenge, and it is remarkable that methods have been developed that allow such things to be done from memory.

Before getting to the code demonstrating different aspects of bell ringing in general and Wolf Wrap in particular, however, we must evaluate the following cell to run the necessary initialization code:

✕

Now, consider the method shown in the following grid:

✕

The three rules are obeyed and the last row is 13527486. Note how the treble bell (#1) moves up from position 1 to position 8 and then back to position 1. The other bells do the same (except for the last change) but with an offset. This pattern is easy for the ringers to remember.

This is a permutation of order 7—the 7-cycle (2, 3, 5, 7, 8, 6, 4):

✕

In the first change in the preceding grid, no bell makes places; that is denoted by ×. At the second change, places are made at positions 1 (by bell 2) and 8 (by bell 7). So the place notation for the sequence above is ×18×18×18×18×18×18×18×12. A ringer would learn this method by memorizing, for each change, the positions at which bells make places (i.e. the bell does not change its position). The previous grid represents the first lead of plain Bob major; “major” indicates eight bells. We will use {} for the ×; thus:

✕

✕

Here is what plain Bob major sounds like. Note the pleasing musical pattern at row 8, which is the reverse of the identity, called back rounds. There is a small pause after every even-indexed row, which is standard practice:

✕

If the place notation is simply repeated seven times, then the final row will be rounds. That is called a plain course of plain Bob major. Here the plain course has seven leads. The permutation that starts each lead is called a lead-head. In the earlier grid, row 16 is the second lead-head; the second lead-head will play an important role in the constraints in §3. There are 7 ⋅ 16, or 112, changes and the plain course sounds like this:

✕

Ringing a peal or quarter-peal is a complex physical and mental exercise: no written music is allowed, so it is possible only because, over several centuries, the ringers developed a concise system for memorizing the patterns that guide the bells. The ringers learn the place notation for the method and the conductor (one of the ringers) will call out certain special changes to be rung at the end of some leads. So while a plain course of plain Bob major has only 112 changes, the special changes allow this to be extended through 5,000 changes with no repetitions.

Rounds—the identity permutation 12345678—are pleasing to listen to but can occur only at the start and finish of a plain course. Mark Davies has devised methods that have several occurrences of rounds, but wrapping around the end of one change to the beginning of the next (this is called a wrap). The following two rows illustrate a wrap:

46587123 45678213

Here, 12345678 occurs at the end of the first row and the beginning of the second. An additional condition is that there should be no pause during the wrap. This means that a wrap should always start in an even-indexed row, using the indexing of the previous grid.My idea was to use ILP to develop a method that maximizes the number of wraps. ILP is available in Mathematica through the LinearOptimization function.

2. A Linear Description of Bell Ringing

To use ILP to invent a new ringing method, we need linear equations and inequalities to represent the basic rules. The starting point is to use variables B(x, y, i) with the conditions 0 ≤ B(x, y, i) ≤ 1 and B(x, y, i) ∊ ℤ. Such a B is called a binary variable: it is 0 or 1. The idea is that a value of 1 for B(x, y, i) means that bell i occurs at grid square (x, y), where y is a row number as indicated in the table of §1. We will restrict to eight bells and use k to denote the number of rows. In what follows, x and i will always vary from 1 to 8. First we set up the variables:

✕

Define constraints, forcing each variable to be 0 or 1 and the first row to be rounds. Also enforce the condition that exactly one bell occupies any position and that bells in any row are a permutation of 1 through 8:

✕

To enforce the critical bell rule 3, we eliminate the possibility of any bell moving more than one position. For example, B(3, 10, 1) + B(5, 11, 1) ≤ 1 says that bell 1 cannot move from position 3 in row 10 to position 5 in row 11:

✕

Place notation plays an important role later, so we show how to force a binary variable P(x, y, i) to indicate whether bell i makes places at position (x, y) (i.e. B(x, y, i)= B(x, y + 1, i)= 1). A natural approach is to set P(x, y, i )= B(x, y, i) ⋅ B(x, y + 1, i), but this is illegal because multiplication is not a linear operation. There is, however, a standard technique in linear programming to reduce such a multiplication to four linear inequalities! I’m grateful to ILP expert Rob Pratt for teaching me this. First, here are the variables, including P(x, y), which indicate whether any bell makes places at the location. A simple sum guarantees that P(x, y) is as desired:

✕

Here is the multiplication trick: four inequalities can be used to force P(x, y, i) to be the product of B(x, y, i) and B(x, y + 1, i). The identity is easily checked by running through the four possibilities for the two B values:

This technique works in a similar way for products of more symbols, and that will be used later in the constraint that counts wraps:

✕

Another commonly enforced condition states that a bell cannot stay in its place three times in a row. This is called no long places. It is easy to set this up in linear fashion using the P variable: P(x, y) and P(x, y + 1) should not both be 1. We work modulo k − 1 so that we can compare the last-place list with the first, as that is relevant:

✕

Truth is tricky and will be discussed later. One subterfuge that often leads to truth is to ask that the place notation does not have ×× and also that it never happens that places are made everywhere:

✕

Here is a check. We use LinearOptimization with objective 0; this means it will search for an assignment of values that satisfies all the constraints. A key point is that for an optimization or feasibility problem with integer variables and linear constraints, LinearOptimization uses ILP. By default, it uses the ILP algorithms of the COIN-OR library. One can also (via the Method function) use the newer and state-of-the-art ILP solver called GUROBI, but a GUROBI license—it’s free for academic use or you can sign up for a free trial—is required.

The following search for a legal method takes well under a second with GUROBI, while the default method takes about 50 seconds:

✕

Here is what LinearOptimization found:

✕

There are no repetitions in the rows, so the conNoDoubleCross constraint worked. Here is the place notation for the preceding; dots are used to separate data for rows when there is no ×:

✕

The last permutation has order 6, so the place notation for the 49-row plain course would consist of six repeats of the place notation:

✕

To skip to the chase, §3 will show how ILP can lead to the discovery of the method described by the following sequence of places

which I called Wolf Wrap. The final permutation is 81234567, which is an 8-cycle (so the plain course has 8 ⋅ 12= 96 changes) and Wolf Wrap has an unusual property that is of interest to the bell-ringing community:

✕

✕

In the wider ringing world, methods often have some symmetry properties. We will not discuss those here; see [2] for more on symmetry and how to implement symmetry constraints linearly.

Example: COIN-OR

Simple examples do not require GUROBI. The default method, COIN-OR, is used here on the example from this section. It is more than one hundred times slower:

✕

3. Wrapping It Up, Truthfully

Our goal is to obtain a method whose plain course has many wraps (and satisfies the truth condition), but we also want the method to have aesthetic and musical appeal. The latter is needed so a band will ring the method; having a successful ring of a quarter-peal or peal is a requirement for a method to become officially named and listed.

To set up an ILP search for a wrap-rich method, we must set up variables to count the number of wraps and also enforce truth (no repeated rows) for the method. These are both a little tricky. Some easier additional constraints are the following; they lead to place notation that is nicely simple for the ringers to memorize:

The number of places made at each change should be either 0 or 2.

The total number of places made in the entire method should not be too large.

These have simple linear implementations:

✕

✕

3.1. Counting Wraps

Both truth and wrap counting are properties of the plain course, not just the first lead. So to set these up as equations, I needed to know the permutation that forms the second lead-head. If that permutation is p (having order m) and if F denotes the first lead, then the plain course consists of F, p(F), p^{2}(F), …, p^{m – 1}(F), where each set is just the set of permutations of each row of F using a power of p. So we will specify that permutation. Looking at some compositions by Mark Davies with high wrap counts led me to set the lead-length (the number of changes in the first lead, denoted L) to 12 and to set the second lead-head to 81234567, which is an 8-cycle. So there are 12 rows in each lead and eight leads in the plain course. Specifying a lead-head permutation is simple:

✕

Recall that in the constraint for places, we made use of a subtle product construction: one can set one binary variable to be the product of two other binary variables via four inequalities. This generalizes to M binary variables as opposed to just two by using 2M inequalities. To get the product W of V_{1}, V_{2}, …, V_{M}, use 0 ≤ W, W ≤ V_{i}, W ≥ 1 − M + ΣV_{i}.

So we can introduce W(x, y, p) to specify that a wrap starts at position (x, y + jL), where x runs from 2 to 8 (1 is excluded because rounds can appear only once), y runs from 1 to 11 in steps of two (because of the parity issue) and p is one of the eight powers of the lead-head permutation. Using the powers in this way provides a method for dealing with a property that can occur in any of the eight leads. The wrap constraint is then the product

W(x, y, p)= B(x, y, p_{1}) ⋅…⋅ B(8, y, p_{8 – x + 1}) ⋅ B(1, y + 1, p_{8 − x + 2})… ⋅ B(x − 1, y + 1, p_{8}),

where the product trick is used to express this using 2 ⋅ 8= 16 inequalities. If this product is 1 and p is the j^{th} power (j= 0, …, 7), then the specified locations have bells 1 through 8 in order, and there is a wrap starting in the (y + jL)^{th} row, starting at position x. The full code for the wrapping and truth constraints is in the supplement.

Recalling the parity issue concerning wrap, the following constraint will then force the number of wraps to equal the target where the initial row of rounds counts as a wrap. This assumes that the lead length is even; it is a little different in the odd case:

✕

3.2. Enforcing Truth

Truth will be handled in a way similar to the wrap count. We will specify the lead-head and then say that the j^{th} row in the first lead is not equal to the m^{th} row in another lead (m ≠ j), using powers of the lead-head as was done for wraps. First, consider how to compare two rows in the first lead. B(x, y_{1}, i) + B(x, y_{2}, i) − 2B(x, y_{1}, i) ⋅ B(x, y_{2}, i)= 1 means that the x position in the two rows are different as regards bell i; a value of 0 means they are the same. So we can simply say that the sum of B(x, y_{1}, i) + B(x, y_{2}, i) − 2B(x, y_{1}, i) ⋅ B(x, y_{2}, i) over all choices of x and i must be at least 1. Of course, we use the product trick discussed earlier (use a new binary variable for each possible product) to linearize the products that arise. And to deal with rows in different leads, we use the permutation that forms the lead-head, which we take as known for our application.

3.3. Discovering Wolf Wrap

So here is the final code to find a method having nine wraps. The lines starting with If[EvenQ[leadLength] where target–1 appears account for the fact that the occurrence of rounds in the first row is considered to be a wrap (it is a legal occurrence of rounds). That line uses the parity of the lead length to choose the y values where a wrap is legal. The parameters here—nine for the target number of wraps and 14 for the total number of places—were obtained after some trial and error led to the conclusion that these values were optimal:

✕

The resulting linear system is enormous, with over 36,000 variables and over 151,000 constraints. It is remarkable that GUROBI finds a solution in about three minutes:

✕

Here is what this method, which I named Wolf Wrap, looks like:

✕

✕

Here is the plain course for Wolf Wrap with the wraps shown in bold, red numbers:

✕

And here is what it sounds like:

✕

The plain course has 96 rows. To turn this into a quarter-peal, switches at various lead-ends must be introduced. Such a composition was constructed by Janet Archibald, who spliced Wolf Wrap with the Stedman triples method. Her composition of a quarter-peal was then rung on June 3, 2021, in just under an hour; see the Ringing Room and Composition Library.

There were 1,272 changes, so the average ringing time for each bell was was 56 ⋅ 60 ∕ (1272 ⋅ ((8 + ½)), or 0.31 seconds. (The ½ is because of the pause that occurs after every two changes.) I am grateful that Mark Davies, who liked the method, asked Janet Archibald to turn it into a quarter-peal composition, and both were members of the band that rang the quarter-peal.

4. ILP Is Awesome

The diversity of problems that can be attacked with ILP is immense because binary variables can be applied to so many optimization problems. I use it in a consulting business where incoming college students present their preferred courses in order, and we assign them to classes so as to optimize the choices from the students’ perspectives. ILP is also used to make up the schedules for professional sports leagues, and there are myriad other applications. Mathematica makes it easy to set up an ILP, and we can now add bell ringing to the list of endeavors where ILP can play a role.

Why are we releasing Version 13 of Mathematica?
As a math and science company, Wolfram doesn’t fear two-digit number integers, especially that particular one between 12 and 14. That is to say, we don’t have triskaidekaphobia, not on Friday the 13th or any other day of the year.
Many nice things can be said about 13 beyond it being a ternary repdigit 1113 or the answer to Prime[6]. It s one of just three known Wilson primes because = 2 834 329, an integer, and is also the fifth Mersenne prime exponent:
#10005
But 13 has much more going for it than just being a prime number, and here are 13 reasons why.
1. First, 13 is part of many linear recurrences. For example, 13 is a Fibonacci number where the limit of ratios between values goes to the golden ratio, phi, 2= + 1:
#10005
#10005
The fifth substitution in the Ammann chair tiling thus has 13 chairs:
#10005
Thirteen is also a tribonacci number, a term coined by Mark Feinberg at the age of 14. The ratio between terms goes to the tribonacci constant, t 3= t 2 + t + 1:
#10005
Tribonacci numbers can be thought of as binary sequences without 000:
#10005
2. Thirteen is also a Narayana cow number, based on psi, 3= 2 + 1, the supergolden ratio:
#10005
The sixth substitution in the psi-quad tiling thus has 13 quads:
#10005
3. There are 13 Archimedean solids, in which both the golden ratio and tribonacci constant feature prominently. Similarly, there are 13 Archimedean duals, all of which make fair dice:
#10005
4. As a number of form n 2 + (n + 1)2, 13 is a centered square number:
#10005
5. A combination lock with three buttons has 13 combinations where every button is pressed once, making it the third Fubini, or ordered Bell, number:
#10005
6. Balanced ternary notation allows three weights and a two-pan scale to measure weights 1 to 13:
#10005
7. The remainder graph of 13 has the following form. To use it, pick a number like 2,522. Then start at 0 and follow 2 blue, 1 red, 5 blue, 1 red, 2 blue, 1 red, 2 blue to find the remainder of dividing 2,522 by 13. This works for any positive integer:
#10005
8. Using steps of N, E, NE or (1,0), (0,1), (1,1), there are 13 ways to get from (0,0) to (2,2), part of a combinatorial sequence known as the central Delannoy numbers:
#10005
These values can be calculated with Legendre polynomials:
#10005
9. There are 13 books in Euclid’s Elements . For example, Book I, Proposition 1 discusses the making of an equilateral triangle:
#10005
10. A sparse ruler of length n allows all integer distances up to n to be measured with a minimal number of marks. The sparse ruler for length 211 uses 26 marks and features 13 spacings of 13:
#10005
11. The BioSequence for the string "THIRTEEN" is as follows:
#10005
12. Books like An Elementary Introduction to the Wolfram Language use ISBN-13 numbers:
#10005
13. Finally, for obscuring text in a 26-letter alphabet, the standard is the ROT13 Caesar cipher that rotates each letter forward 13 places:
#10005
This allows us to finish this article with an obscure piece of code:
#10005
Engage with the code in this post by downloading the Wolfram Notebook
?Комментарии (2)