Here is a tentative answer to your questions. I have an open bounty for a more thorough answer to question $1$), as I find my answer slightly lacking.
First, I want to focus on two key, related elements of the BKT transition: the universal jump in the spin-stiffness/superfluid density of $\frac{2}{\pi}$ and the universal critical exponent $\eta(T_{c}) = \frac{1}{4}$. Here, I define $\eta$ in terms of the decay of the spin-spin correlation function $C(r) \sim \frac{1}{r^{\eta}}$, or equivalently as the dependence of the susceptibility on the correlation length, $\chi = \xi^{2-\eta}$. Note that for the definition of $\eta$ via $C(r)$, I believe $\eta$ is only universal at the critical point $T=T_c$ and otherwise varies with temperature below $T_c$.
Given the universality above, we could imagine a "BKT-like" transition in higher dimensions. If we imagine $\xi \sim e^{\frac{c}{\sqrt{T-T_c}}}$ as usual for BKT, and assume the "singular" part of the free energy goes as $f \sim \xi^{d}$, then it's clear $f$ is infinitely differentiable but nonanalytic regardless of the dimension $d$. The transition will not be BKT but just "BKT-like" if we can further demonstrate that the universal features discussed above (the universal jump and $\eta(T_c)$) differ from those of the BKT class.
The paper "Topological Phase Transitions in Four Dimensions" gives exactly such a "BKT-like" transition in four dimensions. Here, the universal jump appears to be $\frac{4}{\pi^2}$ and $\eta(T_c) = \frac{1}{32}$.
A couple facts about the model: the effective Hamiltonian density in the action appears to have the four-derivative term $(\nabla^2 \phi)^2$ rather than the more common $(\nabla \phi)^2$. This action gives rise to a lower critical dimension of $4$ rather than the usual lower critical dimension of $2$.
I might be misinterpreting their discussion around equations $1$ and $2$, but I take it to mean that this transition is somewhat similar to a tricritical version of BKT. That is, they argue that (at least in mean field theory), their model can be realized as an XY model with nearest-neighbor and next nearest-neighbor interactions, where the ratio of the interaction strengths is appropriately tuned. The extra fine-tuning is needed because of the necessity to cancel out quadratic derivative terms.
This model thus answers your questions 1) and 2) in the affirmative: there are classical infinite order phase transitions and classical continuous phase transitions without symmetry breaking which are not in the BKT universality class. However, I find my answer somewhat unsatisfactory - as I note above, the transition I'm considering is akin to a tricritical BKT. It would be nice if there were an example without extra fine-tuning!
For your question 3), I will quote from the following March Meeting abstract. If someone can find the talk or an associated paper, I would be much obliged.
The key quote is
By considering the infinite p limit of a free energy that we have
derived for a p-th order phase transition, we can derive a Landau type
free energy. We will discuss the properties of the free energy and
identify the features essential for a description of an infinite order
phase transition. These include a logarithmic interaction between the
fields and a novel dependence on spatial gradients. Contrary to
popular belief, since some symmetry is broken at each finite p order,
we submit that an infinite order phase transition does not exclude a
symmetry being broken. Restricting to one dimension, we solve for
domain wall solutions.
That is, this talk appears to say the following: a logarithmic interaction is necessary for infinite order transitions, lack of symmetry breaking is not necessary for infinite order transitions, and 2$d$ is not necessary.
We have more points of evidence for these three statements. 2$d$ isn't necessary, given the 4$d$ example above. The 4$d$ example above indeed has a logarithmic interaction between vortices (see equations B$13$, B$14$, and B$17$ in appendix B of Topological Phase Transitions in Four Dimensions), which is supporting evidence for the importance of such interactions in infinite-order transitions regardless of dimension.
For the case of symmetry-breaking: There are two BKT (or, at the very least, BKT-type, as I'm unsure of the universality of $\eta$) transitions in planar $q$-state clock models with $q$ large enough. At $q\to \infty$, the model becomes the $XY$ model without an ordered phase, but for $q < \infty$, the ordered phase exists in these models since they have only discrete symmetries. For $5 \leq q < \infty$, the phase diagram looks like, from low temperature to high temperature, "ordered phase, critical phase, disordered phase", with two BKT-type transitions separating the middle critical phase from the surrounding ordered and disordered phases. I think it's fine to view the lower-$T$ BKT transition as an ordering transition, since it separates an ordered phase and a critical phase. For a thorough numerical investigation, I recommend the paper "Critical properties of the two-dimensional q-state clock model." Thus, it appears that infinite-order transitions can coincide with symmetry-breaking transitions.
Let's review my answer to question $1$ one more time. Throughout my answer above, I've really given several example models with classical infinite order phase transitions. There's the usual 2$d$ BKT transition, there's the model I'm calling the 4$d$ tricritical BKT-type transition, and there's the BKT-type transition coinciding with ordering in the $q$-state clock model (the BKT-type transition at lower temperature). All three are infinite-order transitions, but they have different $\eta(T_c)$ values, signaling different universality classes. However, it's not so clear to me that these are really truly fundamentally different from BKT. Indeed, the references above really just call all three of these transitions plain "BKT" instead of my moniker of "BKT-like" for the latter two. I have an open bounty for question 1 - I hope there might be other examples given that are perhaps more radically different from BKT!
EDIT: I will note another infinite-order classical phase transition without symmetry breaking that is not in the BKT universality class. This transition is described in a paper investigating the effect of adding phase disorder to the $XY$-model. At low temperatures, there's a transition induced by changing the disorder strength $\sigma$. Here, the correlation length diverges as $\xi \sim e^{\frac{c}{|\sigma - \sigma_c|}}$, and $\eta = \frac{1}{16}$. The more dramatic dependence of the correlation length on the tuning parameter shows that this is in a universality class distinct from BKT. I note that while the transition here is tuned with $\sigma$, it should also be able to be induced by tuning $T$.
