Jelly, 67 bytes

^{Probably way longer than possible so I might try implementing a perimeter-walk approach.}

ạỊẠɗƇ@ạþ`§Ị§3_RƇS<3
,`WW+ⱮØ.ṚNƭƬ¤ẎçƇ⁹ḟɗṭ€Ʋ€Ẏ$’}¡µŒṬZẸƇ$⁺ZṚ$ƬUƬẎṀ)QL

A monadic Link that accepts the number of squares, \$n\ge 1\$, and yields the number of free polyominoes without holes.

Try it online! (Too slow to complete the first case with holes, \$n=7\$).

To show the hole avoidance, here is augmented code that prints the next polyominoes the code rejects, while this prints the next polyominoes the code allows. (Note if trying out other inputs: only positions that are orthogonally adjacent to an existing square are considered for new squares - i.e. this does not reject anything since all orthogonal continuations are valid, specifically row 2 column 3 isn't under consideration.)

How?

ạỊẠɗƇ@ạþ`§Ị§3_RƇS<3 - Helper Link, Won't make a hole?: PotentialCoord; UsedCoords
   ɗƇ@              - filter keep {UsedCoords} for which:
ạ                   -   {UsedCoord} absolute difference {PotentialCoord} (vectorises)
 Ị                  -   insignificant? (vectorises)
  Ạ                 -   all?
                      -> N = UsedCoords in the neighbourhood of PotentialCoord
                             (the "neighbourhood" being the closest 9 coords, including itself)
      ạþ`           - {N} table of absolute difference (vectorises) with {N}
         §          - sums -> a table of all of the Manhatten distances between coords in N
          Ị         - insignificant? (vectorises) -> table of "is immediate neighbour"?
           §        - sums -> counts of immediate neighbours in N for each coord in N
            3_      - subtract from three (vectorises)
              RƇ    - keep those that are positive (keep if range(x) non-empty)
                S   - sum
                 <3 - less than three?

,`WW+ⱮØ.ṚNƭƬ¤ẎçƇ⁹ḟɗṭ€Ʋ€Ẏ$’}¡ - Main Link (part 1): Area
,`WW                         - pair {Area} with itself and wrap twice
                               -> FoundCoordLists = [[[Area, Area]]]
                               (offsets the starting coord to save a little in part2)
                        $’}¡ - repeat this Area - 1 times:
                     Ʋ€      -   for each of {FoundCoordLists}:
            ¤                -     nilad followed by links as a nilad:
      Ø.                     -       literal [0,1]
           Ƭ                 -       collect up while distinct under:
          ƭ                  -         alternate between calling:
        Ṛ                    -           reverse
         N                   -           negate
                                   -> [[0,1],[1,0],[-1,0],[0,-1]]
    +Ɱ                       -     map across those with: {FoundCoordList} add {it}
             Ẏ               -     tighten -> PotentialCoords (includes used ones & repeats)
                  ɗ          -     last three links as a dyad - f(PotentialCoords, FoundCoordsList):
               Ƈ             -       keep those {PotentialCoords} for which:
              ç ⁹            -         call the HelperLink
                 ḟ           -       discard those in {FoundCoordsList}
                   ṭ€        -     tack each onto a copy of FoundCoordsList
                       Ẏ     -   tighten -> new FoundCoordLists

µŒṬZẸƇ$⁺ZṚ$ƬUƬẎṀ)QL - MainLink (part 2) 
µ               )   - for each of {FoundCoordsLists}:
 ŒṬ                 -   convert to a 2d array with 1s at given coords
      $⁺            -   do this twice (to remove empty rows/columns):
   Z                -     transpose
    ẸƇ              -     keep if all
          $Ƭ        -   collect while distinct under (getting all distinct rotations):
        Z           -     transpose
         Ṛ          -     reverse
             Ƭ      -   collect while distinct under (turning those all over):
            U       -     reverse each
              Ẏ     -   tighten -> all versions under the symmetries
               Ṁ    -   maximum -> a canonical form of the FoundCoordsList
                 Q  - deduplicate
                  L - length

Ruby, 221 220 bytes

->n{c=0
(9**-~n).times{|i|r=i.to_s(3).bytes
(0...2*m=r.size).map{|j|[r=r.rotate,r.reverse][j/m]}.max==r&&(v=1
a=x=y=u=0
q=r.map{|k|k!=49&&(u,v=-v*k-=49,u*k);a+=(v*x+=u)-u*y+=v;[x,y]}
x|y|u|~-v|a-n*2==0&&q&q==q&&c+=1)}
c}

Try it online!

Alternate to my previous answer, based on Xnor's perimeter approach, but it works slightly differently. All walks are first generated as a base 3 number in string form, which is then converted into an array of ascii codes. These are rotated and reflected, and discarded if they are not in canonical (lexically highest) form. Hence there is no need to store the polyominoes.

Once a walk in canonical form is found it is checked to see if it is a valid polyomino. The ascii code 49 (for 1) is considered straight and 48 and 50 are (for 0 and 2) considered left and right.

The algorithm for finding the polyomino area is based on a classic way of finding the area of a polygon: select a point (such as the origin) which can be inside the polygon (but doesn't have to be if areas are signed.) Draw a triangle from each edge to the point and the area contribution to the polygon for each edge is found from the area of a triangle formula base*height/2. The code actually calculates double the area and subtracts it from n*2 instead of dividing.

The following expression confirms that we have returned to the origin, facing the same direction, with area equal to n and without visiting the same point twice: x|y|u|~-v|a-n*2==0&&q&q==q. We use the OR operator to confirm multiple expressions are zero with only a single ==0, then we deduplicate list of points q by self intersection before comparing it with q itself to confirm there are no duplicates.

Canonical path detection removes reflections of the path and starts from different points, but certain reversals of direction are not removed - it considers RRSRRS and SRRSRR to be the same but LLSLLS to be different. This does not matter, as one is considered to have positive area and the other negative.

Maple, 283 bytes

proc(n)t:=[1,I,-I];p:=();for m to 2*n+2 do W:=Iterator:-Necklace(m,3);for w in W do z:=seq(w);z=sort([[z],seq([seq([z,z][-i-j],i=1..m)],j=1..m)])[1];c:=0;a:=0;d:=1;C:={for i to m do c+=d;a+=Re(d)*Im(c);d*=t[z[i]+1];c od};if`and`(c=0,d=1,nops(C)=m,a=n)then p,=z fi od od;nops({p})end;

proc(n)
 t:=[1,I,-I];                 # multiplier to turn by 0 degrees (S), 90 (L), -90 (R) in complex plane
 p:=():                       # initialize set of polyominos
 for m to 2*n+2 do            # perimeter lengths (could be 4 to 2*n+2 by 2 but we're golfing)
   W:=Iterator:-Necklace(m,3);# iterator for necklaces of length m with 3 bead types (turns: S,L,R)
   for w in W do              # for each necklace (beads coded with 0..2 for S, L, R)
     z:=seq(w);
     # next line finds lex lowest canonical ordering of reversed and rotated lists
     z:=sort([[z],seq([seq([z,z][-i-j],i=1..m)],j=1..m)])[1];
     c:=0;                    # initialize coordinates (points in complex plane) to zero
     a:=0;                    # initialize area to zero
     d:=1;                    # arbitrary start direction (+x)
     C:={                     # find coordinates to detect self-intersection and current y value
         for i to m do        # for each bead                  
           c+=d;              # update coordinate                                     
           a+=Re(d)*Im(c);    # update area          
           d*=t[z[i]+1];      # decode and update direction
           c                  # last value in loop retained in set C
         od   
        };
     if `and`(c=0,d=1,        # valid if returns to origin in right direction             
            nops(C)=m,        # and no duplicate coordinates (not self intersecting)
            a=n)              # and correct area
     then
       p,=z                   # record this polyomino
     fi          
   od
 od;
 nops({p})                    # remove duplicates and count
end;

This is @xnor's nice algorithm with @Level River St's canonical form idea. Coordinates are coded as complex numbers. Areas are done by the sum of (Delta x)*y formula. Maple has a built-in iterator for necklaces, so only rotationally unique perimeters are generated (in lowest lex order canonical form). Unfortunately, it doesn't have a bracelet iterator, so the reversed lists have to be done separately.

n = 5 takes 7 s, n = 6 takes 76 s.

Python 3, 221 bytes

def f(n):
 L=[[]];c=0
 for l in L:
  A=p=0;s,*r=1,;L+=l+[0],l+[1],l+[3];L*=len(l)<2*n+3
  for d in l:s*=1j**d*2**(p in r);r+=p,;p+=s;A+=p/s
  c+=s-1==p==0<A.imag==n*2!=l==max(max(l:=l[1:]+[x],l[::-1])for x in l)
 return c

Try it online!

Saved 8 bytes with Level River St's idea of discarding paths not in canonical form

I haven't golfed this too hard -- I just want to show off a new method.

We represent a hole-free polyomino by its boundary loop walk. This means, tracing its perimeter one segment at a time, which direction do we turn at each step: left, right, or straight?

For example, the area-3 polyomino

XX
X

has the boundary walk

L<S<L
v   ^
S R>L
v ^
L>L

which corresponds to RLLSLSLL if we start from the R.

This representation has some nice features compared to the usual of a polyomino as the coordinates of each of its squares. Because everything is relative, there's no need to handle translated or rotated polyominoes. Ones with interior holes are never generated, and ones with self-touches can be culled by checking if the boundary ever revisits a point.

Moreover, there's a neat "discrete calculus" trick to find the area from the walk, and I've even posted a challenge about it before: Area enclosed by perimeter loop!

So, here's the method. Generate all sequences of turns (left, right, straight) of sufficiently large length. For each one, check that it:

• Returns to its start point
• Ends facing the starting direction
• Never encounters a coordinate a second time, since this would indicate a loop of self-crossing
• Has the desired area

If so, add it to the list of polyominos. However, polyominos will repeat because the starting point of the loop is arbitrary, and reflections correspond to reversed walks. So, convert each walk to a canonical form, by taking the list-maximum among its cyclic rotations and reversals, and only increment the count if it was already in that form.

This code is slow -- it only gets up to n=5 on TIO. I've tested a somewhat-optimized version up to n=7.

Charcoal, 114 bytes

≔Ｉ1jδＮθ⊞υ⟦⟧Ｆυ«≔ＸδＥιΣ…ι⊕λη≔θζ≔⁰εＦη≡κδ≧⁻εζ±δ≧⁺εζ≧⁺κε¿⁼⁰ΣηＭ∧¬ζ∧⁼¹⊟η⁼ι⌊ΣＥＥιＥι§ι⁺λν⟦⮌κκ⟧→¿›‹Ｌι⊗⊕θ№ＥηΣ…ηλΣηＦ³⊞υ⁺ι⟦⊖κ⟧»Ｉⅈ

Attempt This Online! Link is to verbose version of code. Explanation: Uses @xnor's approach, but even slower so can only manage up to n=4 on ATO.

≔Ｉ1jδ

Get i, the square root of minus one, into a variable, saving a byte.

Ｎθ

Input n.

⊞υ⟦⟧Ｆυ«

Start a breath-first search of walks.

≔ＸδＥιΣ…ι⊕λη

Compute the cumulative sums of turns of this walk in terms of powers of i.

≔θζ

Start calculating the area with n squares left to account for.

≔⁰ε

Start at the origin. (I'm not sure it is necessary to reset the value each time in which case three bytes could be saved by using one of the predefined variables with a numeric value instead.)

Ｆη≡κ

Loop over the direction vectors.

δ≧⁻εζ

If this is a step up then subtract the current x-position from the area left to account for.

±δ≧⁺εζ

If this is a step down then add the current x-position to the area left to account for.

≧⁺κε

Otherwise adjust the current x-position.

¿⁼⁰Ση

See whether the sum of the direction vectors is 0, indicating that this is a closed loop. Conveniently in Charcoal, Sum([]) is None, not 0, so this fails for the the empty walk.

Ｍ∧¬ζ∧⁼¹⊟η⁼ι⌊ΣＥＥιＥι§ι⁺λν⟦⮌κκ⟧→

If this closed loop has the correct area, ends facing the original direction, and is in canonical form, the increment the number of polyominoes found.

¿›‹Ｌι⊗⊕θ№ＥηΣ…ηλΣη

Otherwise if this open loop is short enough and does not self-intersect, then...

Ｆ³⊞υ⁺ι⟦⊖κ⟧

... add left, straight and right to the walk and save the resulting walks to the list of walks to process.

»Ｉⅈ

Output the final number of polyominoes found.

Ruby, 275 bytes

f=->n,p=[[2]],s=1{q=[]
z=[w=n*3,-1,-w,1]
p.map{|a|(s*4).times{|i|r=[];g=h=1
c=a+[d=z[i%4]+a[i/4]]
8.times{|j|h+=(g-g=2*a.count(d-z[j/2]+j%2*z[j/2-1]))**2
c=c.map!{|k|j%2<1?k-k%w*2:k%w*w+k/w}.sort!.map{|m|m-c[0]+n}
r<<c*1}
a.any?(d)||h>10||q<<r.min}}
s<n ?f[n,q&q,s+1]:p.size}

Try it online!

Adapted from my answer to a previous question about polyiamonds (shapes composed of equilateral triangles).

This is a recursive function taking required size of polyominos n as an argument, plus two optional arguments: an array p of smaller polyominos and s the number of squares in the smaller polyominos. The initial default value of these corresponds to the seed polyomino of just one square.

First golfed version, runs a little slower than the first working version below due to certain checks being carried out late in the algorithm to avoid (), resulting in a lot of unnecessary calculation. In particular, all next generation polyominos are calculated before a check is made to see whether we have already reached the required value of n, as the decision to return the correct value or recurse further is made at the end of the function. Similarly detection of new squares overlapping old ones is not done until sufficient info has been gathered to find if a hole has been created.

The polyomino is represented as an array of integers of the form d=x+wy. The board is of width 3*n in the x direction to avoid overlaps and is effectively infinite in the y direction.

The main new feature is the avoidance of holes. The polyomino is grown by adding a new square beside any existing square in any of the four orthogonal directions. The new square is rejected if it bridges across to another existing square. This is detected by walking around the 8 spaces adjacent to the new square and counting the number of changes from full square to empty square and vice versa. There should be exactly 2 changes. If there are more, the new square is bridging and the polyomino is rejected.

The golfed version maintains a value of 0 or 2 in g for the previous square and calculates a new value each iteration, adding the square of the difference to h each time. h therefore increments by 4 each time a change is detected. Because there are an even number of changes the idealised value of h would be a multiple of 8. In the golfed version we take advantage of this by initialising g and h to a fuzzy value of 1. The first increment of h is therefore always to 2, followed by 1-4 proper change detections. If there are more than 2 proper change detections, h exceeds 10 and the polyomino is rejected.

At the same time, the polyomino is alternately reflected 8 times orthogonally (by subtracting x*2=k%w*2) and diagonally (by swapping x and y k%w*w+k/w) to generate all 8 possible orientations.

If the polyomino is not rejected, it is added to the new list of larger polyominos q.

Commented code

f=->n,p=[[2]],s=1{                  #Call function with desired n, seeded with single square
q=[]                                #Make an array for next generation polyominos
z=[w=n*3,-1,-w,1]                   #w=width of board. z=array of 4 orthogonal directions
p.map{|a|(s*4).times{|i|            #Try extending each existing polyomino in p, from each square in 4 directions
  r=[];g=h=1                         #Make array to hold all orientations of new shape. Initialise g&h for hole detection
  c=a+[d=z[i%4]+a[i/4]]              #Make new polyomino c by extending from a[i/4] in direction i%4 adding square d  
  8.times{|j|h+=                          #Loop 8 times for relections and hole detection
(g-g=2*a.count(d-z[j/2]+j%2*z[j/2-1]))**2 #Walk around 8 neighbours of square d, update g, and if it changed increment h
    c=c.map!{|k|j%2<1?k-k%w*2:k%w*w+k/w}. #Reflect c orthogonally or diagonally depending on value of j%2
    sort!.map{|m|m-c[0]+n}                #sort & normalise the elements so that the lowest=n putting it in middle of board
    r<<c*1}                               #add current relection of c to r. We will select the lexically lowest value of c
  a.any?(d)||h>10||q<<r.min}}        #If new square d doesn't overlap existing square & h shows no bridge, add new shape to q
s<n ?f[n,q&q,s+1]:p.size}           #if s==n return the size of existing list p. 
                                    #If s<n, take precalculated list of next size polyominos q, eliminate duplicates and recurse

Ruby, 350 bytes

f=->n,p=[[2]],s=1{z=[w=n*3,w-1,-1,~w,-w,1-w,1,-~w];p p.size
#p.map{|a|t=(?.*~-w+$/)*s;a.map{|k|t[k]=?#};puts t,$/}
s<n&&(q=[]
p.map{|a|(s*4).times{|i|d=z[i%4*2]+a[i/4]
a.index(d)||(c=a+[d];r=[];g=a.count(d-~w);h=0
8.times{|j|c=c.map!{|k|j%2<1?k-k%w*2:k%w*w+k/w}.sort!.map{|m|m-(c[0])+n};r<<c*1
h+=g^g=a.count(d+z[j])}
h<3&&q<<r.min)}}
f[n,q&q,s+1])
}

Try it online!

First working version This version prints the full sequence from 1 to n. It will even print the polyominos themselves if you uncomment the second line.

Python, 686 bytes

Python port of formatted Ruby code by @Level River St

Golfed version. Attempt This Online!

H=range
def G(n,p=None,s=1):
    if p is None:p=[[2]]
    A=n*3;I=[A,A-1,-1,~A,-A,1-A,1,-~A];print(len(p))
    if s<n:
        C=[]
        for D in p:
            for J in H(s*4):
                E=I[J%4*2]+D[J//4]
                if E not in D:
                    B=D+[E];K=[];F=D.count(E-~A);L=0
                    for M in H(8):
                        if M%2<1:B=[B-B%A*2 for B in B]
                        else:B=[B%A*A+B//A for B in B]
                        B.sort();B=[A-B[0]+n for A in B];K.append(B.copy());N=F;F=D.count(E+I[M]);L+=N^F
                    if L<3:O=min(K);C.append(tuple(O))
        C=list(set(C));C=[list(A)for A in C];return G(n,C,s+1)
    else:return p

Ungolfed version. Attempt This Online!

def f(n, p=None, s=1):
    if p is None:
        p = [[2]]
    w = n * 3
    z = [w, w - 1, -1, ~w, -w, 1 - w, 1, -~w]

    print(len(p))  # Equivalent to `p p.size` in Ruby

    if s < n:
        q = []
        for a in p:
            for i in range(s * 4):
                d = z[(i % 4) * 2] + a[i // 4]
                if d not in a:
                    c = a + [d]
                    r = []
                    g = a.count(d - ~w)
                    h = 0
                    for j in range(8):
                        # Transform `c` based on the value of `j`
                        if j % 2 < 1:
                            c = [k - (k % w) * 2 for k in c]
                        else:
                            c = [(k % w) * w + (k // w) for k in c]

                        c.sort()
                        c = [m - c[0] + n for m in c]
                        r.append(c.copy())

                        # Update `h` with XOR of previous `g` and new `g`
                        temp = g
                        g = a.count(d + z[j])
                        h += temp ^ g

                    if h < 3:
                        min_r = min(r)
                        q.append(tuple(min_r))  # Use tuple to make it hashable for set operations

        # Remove duplicates by converting to a set and back to a list
        q = list(set(q))
        q = [list(item) for item in q]  # Convert tuples back to lists

        return f(n, q, s + 1)
    else:
        return p


# Execute the function with n = 12
result = f(12)