J-uby, 58 bytes

:+&1|:*|:zip%(:*&(~:**&0.5|Z))|:*&+:-|~:each_cons&2|:*&+:+

Attempt This Online!

Japt, 8 bytes

Another port of xnor's formula.

ôÈ-¬fÃä+

Try it

ôÈ-¬fÃä+     :Implicit input of integer
ô            :Range [0,input]
 È-          :Map & subtract
   ¬         :  Square root
    f        :  Floored
     Ã       :End map
      ä+     :Consecutive pairs reduced by addition

Husk, 10 bytes

ΘẊ+mṠ-o⌊√ḣ

Try it online!

Jelly, 8 bytes

R_ƽ$+ƝŻ

Try it online!

-5 bytes by porting xnor's formula (thanks Bubbler!)

-1 byte thanks to Nick Kennedy

Explanation

Uses xnor's formula of:

$$f(n) = \sum_{k \in \{n,n+1\}}\left({k-\lfloor \sqrt k\rfloor}\right)$$

R_ƽ$+ƝŻ    Main Link
R           range
 _   $      subtract
  ƽ                 square root floored (of each element)
     +Ɲ     add adjacent pairs together
       Ż    prepend 0

Without xnor's formula, I have 10 bytes

Jelly, 10 bytes

RƲẸ$Ɲ¬‘ÄŻ

Try it online!

(range; for each pair of adjacent elements, check if either of them is square; logical NOT that and add one (gets the original 1,2 sequence), cumulative sum, prepend 0)

Jelly, 7 bytes

Ż_ƽ$+Ɲ

A monadic Link accepting a positive integer, n, which yields a list of the first n entries.

Try it online!

How?

Application of xnor's pair-wise addition formula \$f(n) = \sum_{k \in \{n,n+1\}}\left({k-\lfloor \sqrt k\rfloor}\right)\$

Ż_ƽ$+Ɲ - integer, n                       e.g. 10
Ż       - zero range                            [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
    $   - last two links as a monad:
  ƽ    -   integer square-root (vectorises)    [0, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3]
 _      -   subtract (vectorises)               [0, 0, 1, 2, 2, 3, 4, 5, 6, 6, 7]
      Ɲ - for neighbours:
     +  -   add                                  [0, 1, 3, 4, 5, 7, 9, 11,12,13]

Ruby, 52 50 bytes

->n,*w{n.times{|x|w+=[x,x*x,x*x];p w.sort[x+1]+x}}

Try it online!

How:

The difference between n and f(n) shows an interesting pattern:

n   f   f-n
----------
0   0   0
1   1   0
2   3   1
3   4   1
4   5   1
5   7   2
6   9   3
7   11  4
8   12  4
9   13  4
10  15  5
11  17  6
12  19  7
13  21  8
14  23  9
15  24  9
16  25  9
17  27  10
18  29  11
19  31  12

In the rightmost sequence, every non-square number occurs only once, and every square number appears three times (except 0 which occurs only twice). I can build the required sequence as the sum of n and (f-n).

C (gcc), 81 bytes

r,c,i;f(n){for(r=c=0;~n;c++)for(i=++c;i--+2&&n--;r+=c+~i&&i+2)printf("%d ",r++);}

Some magic hackery used.

-4 bytes thanks to ceilingcat!

Try it online!

C (gcc) -lm, ^{72 \$\cdots\$ 61} 60 bytes

Saved a byte thanks to ceilingcat!!!

s;i;f(n){for(s=i=0;i<n;)printf("%d ",i-~i-s-(s=sqrt(++i)));}

Try it online!

Uses xnor's formula.

Perl 5 -MList::Util=sum -n, 39 bytes

Shoutout to @xnor for the formula. This is essentially a port of his Python answer.

map{say$a+($n=$_-int$_**.5);$a=$n}1..$_

Try it online!

Charcoal, 17 bytes

ＩＥＥＮ…±⊕ι⊕ιＬ⁻↔ιＸι²

Try it online! Link is to verbose version of code. Based on @xnor's formula. Explanation:

   Ｎ                Input as a number `m`
  Ｅ                 Map over implicit range `0`..`m-1`
       ι ι          Current index `n`
      ⊕ ⊕           Incremented (i.e. `1`..`m`)
     ±              Negated
    …               Exclusive range (i.e. `-n` .. `n-1`)
 Ｅ                  Map over list of ranges
             ι ι    Current range
              Ｘ ²   Squares of values
            ↔       Absolute values
           ⁻        Remove the squares
          Ｌ         Take the length
Ｉ                   Cast to string
                    Implicitly print

JavaScript (ES7), 40 bytes

This is using the closed-form formula described below.

But because we're asked to output the \$n\$ first terms of the sequence, we need 19 bytes of wrapping code. :'-(

f=n=>n?[...f(n-1),(n-=n**.5)*2|n%1>0]:[]

Try it online!

21 bytes (n-th term, 1-indexed)

n=>(n-=n**.5)*2|n%1>0

Try it online!

Given \$n\ge0\$, we compute:

$$d(n)=2\cdot\lfloor n-\sqrt{n}\rfloor\\ f(n)=\cases{ d(n)&\text{if $n$ is a square}\\ d(n)+1&\text{otherwise} }$$

The JS implementation uses a bitwise OR which implicitly floors \$n-\sqrt{n}\$ after it has been multiplied by \$2\$. But this leads to the same result.

05AB1E, 9 7 bytes

ÝDtï-ü+

Port of @Bubbler's top APL answer, which uses the same formula as @xnor's Python answer:
$$f(n) = \sum_{k \in \{n,n+1\}}\left({k-\lfloor \sqrt k\rfloor}\right)$$

-2 bytes thanks to @Grimmy.

Try it online.

Explanation:

Ý        # Push a list in the range [0, (implicit) input-integer]
 D       # Duplicate this list
  t      # Take the square-root of each value
   ï     # Cast it to an integer to floor it
    -    # Subtract the values at the same positions from one another
     ü   # For each overlapping pair:
      +  #  Add them together
         # (after which the result is output implicitly)

Implementing the steps described in the challenge description would be 13 bytes instead:

2∞и1δš€û˜.¥I£

Try it online.

Or 2∞и1δš€û could alternatively be ÅÉÅ21δ.ø.

Try it online.

Explanation:

 ∞             # Push an infinite positive list: [1,2,3,...]
2 и            # Repeat 2 that many times as list: [[2],[2,2],[2,2,2],...]
    δ          # For each inner list:
   1 š         #  Prepend a leading 1: [[1,2],[1,2,2],[1,2,2,2],...]
      €        # For each inner list:
       û       #  Palindromize it: [[1,2,1],[1,2,2,2,1],[1,2,2,2,2,2,1],...]
        ˜      # Flatten the list of 1s and 2s: [1,2,1,1,2,2,2,1,1,2,2,2,2,2,1,...]
         .¥    # Undelta it (cumulative sum with 0 automatically prepended):
               #  [0,1,3,4,5,7,9,11,12,13,15,17,19,21,23,24,25,...]
           I£  # Leave the first input amount of items from this infinite list
               # (after which the result is output implicitly)

ÅÉ             # Push a list of odd numbers below or equal to the (implicit) input
               #  i.e. 6 → [1,3,5]
  Å2           # Repeat a list of 2s for each inner value: [[2],[2,2,2],[2,2,2,2,2]]
     δ         # For each inner list:
    1 .ø       #  Surround it with 1s: [[1,2,1],[1,2,2,2,1],[1,2,2,2,2,2,1]]
               # (The rest is the same as above)

05AB1E, 9 bytes

ENŲ_©O=®

Try it online!

E           # loop for N from 1 to input:
 NŲ        #  is N a square?
    _       #  logical not (0 if N is a square, 1 if not)
     ©      #  save in the register without popping
      O     #  sum all numbers on the stack
       =    #  print without popping
        ®   #  push the content of the register

Python 2, 71 bytes

f=lambda n,k=0,w=3:n*[n]and[0]+[x-~(k>1)for x in f(n-1,~-k%w,w+2*0**k)]

Try it online!

Bash + GNU utilities, 33 29 bytes

seq -f %0.fddv-r1-dv-+p $1|dc

Try it online!

This is another solution using @xnor's nice formula.

Haskell, 42 bytes

(`take`q 4)
q k=0:[1,3..k]++map(k+)(q$k+4)

Try it online!

Uses a version of Bubbler's observation that the sequence alternates runs of consecutive odd numbers with an even number directly in between.

Haskell, 43 bytes

(`take`scanl(+)0(q[2]))
q r=1:r++1:q(2:2:r)

Try it online!

Generates an infinite list of 1's and 2's, take the cumulative sums, and truncates to the input length.

Python 2, 52 bytes

n=p=0
exec"n+=1;r=n-n**.5//1;print p+r;p=r;"*input()

Try it online!

54 bytes

lambda N:[n-~n-n**.5//1-(n+1)**.5//1for n in range(N)]

Try it online!

It's a formula!

$$f(n) = 2n+1 - \lfloor \sqrt n\rfloor - \lfloor \sqrt {n+1} \rfloor$$

This can also be split up as

$$f(n) = \sum_{k \in \{n,n+1\}}\left({k-\lfloor \sqrt k\rfloor}\right)$$

Note that \$k-\lfloor \sqrt k\rfloor\$ is the number of non-squares from \$1\$ to \$k\$ inclusive.

APL (Dyalog Extended), ¹⁴ 12 bytes

0,2+/⍳-⌊∘√∘⍳

Try it online!

Uses xnor's formula of

$$ f(n) = \sum_{k \in \{n,n+1\}}\left({k-\lfloor \sqrt k\rfloor}\right) $$

How it works

0,2+/⍳-⌊∘√∘⍳
     ⍳-       ⍝ 1..n minus...
       ⌊∘√∘⍳  ⍝ floor(sqrt(1..n))
  2+/         ⍝ Add two consecutive pairs
              ⍝ giving first n items of the sequence except leading 0
0,            ⍝ Prepend the leading 0

APL (Dyalog Extended), 14 bytes

⊢↑2(∧+/,2××/)⍳

Try it online!

Based on the observation that the sequence is the union of all odd numbers and the numbers in the form of \$2n(n+1), n \ge 0\$. Uses ⎕IO←0.

How it works

⊢↑2(∧+/,2××/)⍳  ⍝ Input: positive integer n
             ⍳  ⍝ Generate 0..n-1
  2(      ×/)   ⍝ Pairwise product (0×1, 1×2, ..., (n-2)×(n-1))
        2×      ⍝ Double it
     +/,        ⍝ Concat with pairwise sum (0+1, 1+2, ..., (n-2)+(n-1))
    ∧           ⍝ Ascending sort the 2(n-1) numbers in total
⊢↑              ⍝ Take the first n numbers
                ⍝ For n=1, "overtake" from zero elements, giving single 0

Python 3, 93 bytes

f=lambda n,x=0:(n-x)*[1]and[sum([j for i in range(1,n,2)for j in[1]+[2]*i+[1]][:x])]+f(n,x+1)

Try it online!

-21 bytes thanks to @Bubbler