Deranged !Combinatorics: Compute the Subfactorial The subfactorial or rencontres numbers (A000166) are a sequence of numbers similar to the factorial numbers which show up in the combinatorics of permutations. In particular the nth subfactorial !n gives the number of derangements of a set of n elements. A derangement is a permutation in which no element remains in the same position. The subfactorial can be defined via the following recurrence relation:
!n = (n-1) (!(n-1) + !(n-2))

In fact, the same recurrence relation holds for the factorial, but for the subfactorial we start from:
!0 = 1
!1 = 0

(For the factorial we'd have, of course, 1! = 1.)
Your task is to compute !n, given n.
Rules
Like the factorial, the subfactorial grows very quickly. It is fine if your program can only handle inputs n such that !n can be represented by your language's native number type. However, your algorithm must in theory work for arbitrary n. That means, you may assume that integral results and intermediate value can be represented exactly by your language. Note that this excludes the constant e if it is stored or computed with finite precision.
The result needs to be an exact integer (in particular, you cannot approximate the result with scientific notation).
You may write a program or a function and use any of the standard methods of receiving input and providing output.
You may use any programming language, but note that these loopholes are forbidden by default.
This is code-golf, so the shortest valid answer – measured in bytes – wins.
Test Cases
n     !n
0     1
1     0
2     1
3     2
4     9
5     44
6     265
10    1334961
12    176214841
13    2290792932
14    32071101049
20    895014631192902121
21    18795307255050944540
100   34332795984163804765195977526776142032365783805375784983543400282685180793327632432791396429850988990237345920155783984828001486412574060553756854137069878601
Funciton, 336 bytes

Byte count assumes UTF-16 encoding with BOM.

┌─╖┌─╖  ┌─╖ 
│f╟┤♭╟┐┌┤♭╟┐
╘╤╝╘═╝├┘╘═╝├────┐
 │┌─╖ │ ┌┐┌┘╔═╗╓┴╖
 ││f╟─┴┐└┴┼─╢0║║f║
 │╘╤╝  │  │ ╚═╝╙─╜
 │┌┴╖ ┌┴╖┌┴╖ ╔═╗
 ││+╟┐│×╟┤?╟┐║1║
 │╘╤╝│╘╤╝╘╤╝┘╚╤╝
 └─┘ └─┘  └───┘


This defines a function f which takes one integer and outputs another integer at a 90 degree turn to the left. It works for arbitrarily large inputs.

Try it online!

Considering this is Funciton it's even reasonably fast (n = 20 takes about 14 seconds on TIO). The main slowdown comes from the double-recursion, as I don't think the Funciton interpreter automatically memoises functions.

Unfortunately, some monospaced fonts don't correctly monospace the ♭ and/or insert small gaps between the lines. Here's a screenshot of the code from TIO in all its beauty:



I think it might be possible to golf this some more, e.g. by changing the condition from &gt;0 to &lt;1 and swapping the branches of the conditional, so that I could reuse the number literal, or maybe by using a completely different formula, but I'm quite happy with how compact it is already.

Explanation

This basically implements the recursion given in the challenge, although it uses the base case !(-1) = !0 = 1. n-1 and n-2 are computed with the predecessor function ♭, and the intermediate result n-1 is reused in three places. There isn't much more to it, so I'll just quickly go through the control flow:

               ─┐
               ╓┴╖
               ║f║
               ╙─╜


This is the function header which emits the function's input n long the attached line. This immediately reaches T-junction, which simply duplicates the value.

        ┌┐┌┘╔═╗
        └┴┼─╢0║
          │ ╚═╝


The 0 box is just a numeric literal. A 4-way junction computes two functions: the path that leads off the bottom computes 0 &lt; n, which we'll use to determine the base case. The path that goes left separately computes 0 &lt;&lt; n (a left-shift), but we discard this value with the Starkov construct.

         ┌┴╖ ╔═╗
         ┤?╟┐║1║
         ╘╤╝┘╚╤╝
          └───┘


We lead this into the three-way conditional ?. If the value is false, we return the constant result 1. The loose end right of the ? is the function output. I'm twisting it around by 180 degrees here, so that the relative orientation of input and output of f is more convenient in the rest of the program.

If the condition was true, then the other value will be used. Let's look at the path that leads to this branch. (Note that Funciton's evaluation is actually lazy so that this branch will never be evaluated if it's not needed, which makes the recursion possible in the first place.)

        ┌─╖ 
      ┐┌┤♭╟┐
      ├┘╘═╝
      │
     ─┴┐


In the other branch we first compute n-1 and then split the path twice so we get three copies of the value (one for the coefficient of the recurrence, one for the first subfactorial, the last for n-2).

┌─╖┌─╖
│f╟┤♭╟
╘╤╝╘═╝
 │┌─╖
 ││f╟
 │╘╤╝
 │┌┴╖
 ││+╟
 │╘╤╝
 └─┘ 


Like I said, we decrement one copy again with another ♭, then we feed both n-1 and n-2 recursively to f and finally add the two results together in the +.

       ┐
       │
      ┌┴╖
     ┐│×╟
     │╘╤╝
     └─┘


All that's left is to multiply n-1 by !(n-1) + !(n-2).