In different situations,  when there are other functions or a variable appears more than once,  we can use, for example, the software R.  Let's see an example.
A famous manual on the uncertainty of physical measurements (published in 1996 and used for many more years), while in previous situations it proceed by minimums and maximums, in these cases it states that differential calculation techniques must be used.  It gives the following example (referring to a trolley that descends along an inclined plane):  A = L²/(2·S)·(1/T2²−1/T1²)  where  L = 5±0.05, S = 100±0.2, T1 = 0.054±0.001, T2 = 0.031±0.001.  With differential calculation techniques this is obtained:  A = 87.3±7.8.
Let's see what we get with R using EPS(q) that generates random values from −q to q (random number generators have been in use, widely, since 1980; see here):

source("http://macosa.dima.unige.it/r.R")     # I load a library
L = 5+EPS(0.05); S = 100+EPS(0.2); T1 = 0.054+EPS(0.001); T2 = 0.031+EPS(0.001)
A = L^2/(2*S)*(1/T2^2-1/T1^2); range(V)
#  75.96997  99.57375

A is between 76 and 100. Since (100+76)/2 = 88 and (100−76)/2 = 12, I can say that A = 88±12.
This interval is very different from 87.3±7.8!
It is an understandable and more reliable method.

Another example.  T = 3.2±0.05, F(x) = 2^x−x^2, F(T) = ?.

F = function(x) 2^x-x^2
range( F(3.2+EPS(0.05)) )
#  -1.050602  -1.045944

F(T) is between −1.050 and −1.046. Since (−1.050−1.046)/2 = −1.048 and (1.050−1.046)/2= 0.002, I can say that F(T) = −1.048±0.002.
In the case of a single variable function I can also do:

RANGE(F, 3.2-0.05, 3.2+0.05)
# ~   min     F(min)     Max      F(Max)
#  3.212433 -1.050602  3.150000 -1.045944

Obviously (if you have downloaded R) also the cases dealt with the script "interval arithmetic" can be faced up in this way;  x = 0.9±0.05, y in [1.39,1.40]: x = 0.9+EPS(0.05); y = 1.395+EPS(0.005); range(x^y)
#  0.7965033  0.9311825     0.865±0.070