I was unaware from your (very original) post that you required

1. The capability to handle symbols as the range end-points
2. The capability to handle infinities as the range end-points
3. The capability to determine a range from a single inequality

I note with interest that your own posted 'solution' handles none of these - so I guess we can agree that by your own criteria, it is pretty much garbage

I have modified my solution (again) to handle the issue of generating a range  from single inequalities. It seems to wrk for the test case inculded in the code - but if you find a case where it doesn't work, then I respectfully suggest that you fix it yourself.

Can't be that hard - right: and I get annoyed when the OP keeps changing the original problem

f:=proc(expr, var)
        local ulim, llim, cond, j, p;
        if   numelems(expr)=1
        then cond:= `or`( op(0, expr[])=`<=`,
                          op(0, expr[])= `<`
                         ):
             if   op(1, expr[])= var
             then if cond
                  then llim:= -infinity:
                       ulim:= op(2, expr[]):
                  else llim:= op(2, expr[]):
                       ulim:= infinity:
                  fi;
             else if cond
                  then llim:= op(1, expr[]):
                       ulim:= infinity:
                  else llim:= -infinity:
                       ulim:= op(1, expr[]):
                  fi
             fi:
        else for j from 1 by 1 to 2 do
                 cond:= `or`( op([j,0], expr)=`<=`,
                              op([j,0], expr)=`<`
                            ):
                 if   op([j,1],expr)=var
                 then if   cond  
                      then ulim:=op([j,2], expr)
                      else llim:=op([j,2], expr)
                      fi;
                 else if   cond
                      then llim:=op([j,1], expr)
                      else ulim:=op([j,1], expr)
                      fi;
                 fi;
             od:
        fi;
        p:= x-> `if`( type( x, infinity),
                      x,
                      evalf(x)
                    );
        return p(llim)..p(ulim)
   end proc:
conds:=[ { theta <= Pi-arccos(-3/4+(1/4)*sqrt(13+16*sqrt(2))), arccos(3/4) < theta},

         { arccos(-3/4+(1/4)*sqrt(13+16*sqrt(2)))+Pi <= theta, theta < 2*Pi-arccos(3/4)},

         { theta <= 3*Pi-arccos(-3/4+(1/4)*sqrt(13+16*sqrt(2))), 2*Pi+arccos(3/4) < theta},

         { theta >= z, theta< zz},
         
         { theta >= -infinity, theta<Pi/4},

         { theta > -Pi/4, theta < infinity },

         { theta>=Pi},

         { theta < Pi},

         { Pi < theta},

         { Pi > theta}
       ]:
map(f, conds, theta)[];

I forgot to allow for infinities - although it is trivial to fix. My only excuse is that I was trying to avoid the use of max, min etc or indeed anything which relied on 'evaluting' the end points.

The case {theta>=z, theta<zz} *ought* to return z..zz, even if z, zz cannot be evaluated to numbers.

So far as I can tell, no other solution posted handles this case?

Anyhow revised version below which handles the possibility of infinite (as well as unknown) boundaries.

restart:
f:=proc(expr, var)
        local ulim, llim, cond, j, p;
        for j from 1 by 1 to 2 do
            cond:= `or`( op([j,0], expr)=`<=`,
                         op([j,0], expr)=`<`
                       ):
            if   op([j,1],expr)=var
            then if   cond  
                 then ulim:=op([j,2], expr)
                 else llim:=op([j,2], expr)
                 fi;
            else if   cond
                 then llim:=op([j,1], expr)
                 else ulim:=op([j,1], expr)
                 fi;
            fi;
        od:
        p:= x-> `if`( type( x, infinity),
                      x,
                      evalf(x)
                    );
        return p(llim)..p(ulim)
   end proc:


conds:=[ {theta <= Pi-arccos(-3/4+(1/4)*sqrt(13+16*sqrt(2))), arccos(3/4) < theta},

         {arccos(-3/4+(1/4)*sqrt(13+16*sqrt(2)))+Pi <= theta, theta < 2*Pi-arccos(3/4)},

         {theta <= 3*Pi-arccos(-3/4+(1/4)*sqrt(13+16*sqrt(2))), 2*Pi+arccos(3/4) < theta},

         { theta >= z, theta< zz},
         
         { theta >= -infinity, theta<Pi/4},

         { theta > -Pi/4, theta < infinity }
       ]:

map(f, conds, theta)[];

You PDE is second order in 'r' and second order in 'z', whihc means you need four boundary/initial conditions. Since it is mathematically impossible for matlab to provide a solution for such an equation , with three boundary conditions, stop trying to impress. Failing that, post the matlab code which provides the solution, so I can examine.

'regexps' are one of these things which have been around for ages, just not vey well publicised. Whilst I wouldn't want to dissuade from a book purchase, you'd be amazed what would turn up if you just typed  'regexp tutorial' in your friendly search engine

As you have observed, my solution works provided that in the "initial, approximate" solution, one has Digits=10. so perhaps I should have witten this code as

restart;
Digits:=10:
eqn:=2*tan(Pi*t^2)-tan(Pi*t)+tan(Pi*t)*tan(Pi*t^2)^2 = 0:
f:=unapply(lhs(eqn), t):

ans10:=Array(1..14):
a:=-2.01:
for j from 1 to 14 do
       ans10[j]:=RootFinding:-NextZero(f, a);
       a:=ans10[j]:
od:
rr:= seq
      ( [ (ans10[j-1]+ans10[j])/2,
          (ans10[j]+ans10[j+1])/2
        ],
        j=2..numelems(ans10)-1
      ):
rr:= [  [ ans10[1]-(rr[1][2]-ans10[1]), rr[1][2] ],
         rr,
         [ rr[12][2], ans10[14]+ans10[14]-rr[12][2] ]
      ]:
Digits:=50:
seq( fsolve( eqn, t, fulldigits, rr[j][1]..rr[j][2]), j=1..14);

I don't have any real difficulty with the two-stage process

1. Calculate an approximate solution
2. Using the approximate solutions obtained in (1) above, calculate an accurate solution

I merely wished to observe that this process would "work" for the OP's original problem - not that it was necessarily the *best* way to solve the problem

Probably a bit pointless now, but the following works up to Digits=50 - didn't try any further

restart;
eqn:=2*tan(Pi*t^2)-tan(Pi*t)+tan(Pi*t)*tan(Pi*t^2)^2 = 0:
f:=unapply(lhs(eqn), t):

ans10:=Array(1..14):
a:=-2.01:
for j from 1 to 14 do
       ans10[j]:=RootFinding:-NextZero(f, a);
       a:=ans10[j]:
od:
rr:= seq
      ( [ (ans10[j-1]+ans10[j])/2,
          (ans10[j]+ans10[j+1])/2
        ],
        j=2..numelems(ans10)-1
      ):
rr:= [  [ ans10[1]-(rr[1][2]-ans10[1]), rr[1][2] ],
         rr,
         [ rr[12][2], ans10[14]+ans10[14]-rr[12][2] ]
      ]:
Digits:=50:
seq( fsolve( eqn, t, fulldigits, rr[j][1]..rr[j][2]), j=1..14);

I answered this problem once before - see this link

http://www.mapleprimes.com/questions/219396-Multiple-Assignment

where the code I provided executes and provides all the answers you could expect.

You have asked essentially the same question now, several times, and have ignored the answers provided by myself and others. I suggest you check out the answer you received at the above link, because I won't be fixing it again - and I doubt if anyone else will either

Well it seems to work now, but I can;'t help feeling that having to click the 'cancel' button in the 'file upload' pop-up, in order to upload a file (or link) is a little counter-intuitive.

Pretty sure there used to be a 'submit' button - which made a lot more sense

The following *slight* variation on my original code stores a matrix whose first column is t-values, and whose second column is y-yalues - no variable names just numeric values. Once agin you will have to modify the filepath names to somthing appropriate for your installation
 

  restart;
  n := 40:
  h := 40:

  sys1 := [ diff(x(t), t) = 2*t+1,
                diff(y(t), t) = 1+5*t,
                x(0) = 1,
                y(1) = 0
             ]:
#
# Define a vector in which to place solutions
# OP only has three solutions, but I suppose
# might want more - define variable
#
   nsols:=3:
   sols:=Vector(nsols):
   for j from 1 by 1 to nsols do
   #
   # Compute general solution
   #
      res := dsolve(sys1, numeric);
   #
   # Get solution at the specified argument
   #
      sols[j]:=res(h*0.5*(j-1)/(2*n));
   end do:
   matSols:=convert
                  ( [ seq
                      (  [ rhs(sols[j][1]),
                           rhs(sols[j][3])
                         ],
                         j=1..nops(sols)
                      )
                   ],
                   Matrix
                 );
#
# Export solution matrix to a file. NB, OP will
# have to use something appropriate for his/her
# machine
#
   ExportMatrix
   ( "D:/Users/TomLeslie/myMaple/tryMat",
      matSols
    );
#
# Restart, clear memory etc and then import the
# solution which was previously exported, just
# to check
#
   restart;
   ip:=ImportMatrix
        ( "D:/Users/TomLeslie/myMaple/tryMat" );

Your original post stated (my emphasis)

"I have in mind all the real roots of the equation"

I believe that the code I supplied provides all the real roots in the range [-2..2]. My choice of the range for the "real" roots was governed by the expression which you used - as in (quoted from your original)

a := Student[Calculus1]:-Roots(2*tan(Pi*t^2)-tan(Pi*t)+tan(Pi*t)*tan(Pi*t^2)^2 = 0, t = -2 .. 2)

The code I provided supplies all the real roots in the range [-2..2] - that is what you asked for, and that is what you got

I'm pretty sure this is what Preben has already suggested: stages are as follows

1. Solve your original BVP between [-0.5, 0]
2. From this solution, generate values for f(0),D(f)(0), D@@2(f)(0), theta(0), and D(theta)(0). You now have five 'initial' values which you can use to solve the same problem from 0 to any positive value
3. Similarly from the solution to the original BVP solution at (1) above, generate values for f(-0.5),D(f)(-0.5), D@@2(f)(-0.5), theta(-0.5), and D(theta)(-0.5). You can now use these as 'initial' values to solve the IVP from -0.5 to any, more negative value.
4. Stitch all three solutions together on the same graph.
5. I would have tried this, but since you have not posted your code in any usable way (and I'm not going to retype it!!) - I can't

Do I like this approach? - no, not really: IVPs are often very sensitive to their initial values, so using sets of 'approximate initial values' from the output of a BVP problem is inherently risky. If the solution for the original BVP is reasonably 'smooth' and 'bounded' on the boundaries, then this approach might be viable.

Certainly can't think of any other way to do it!!!

My only excuse is that I was fighting the update to the MaplePrimes interface, and could no longer insert worksheets; so I just typed without checking - my bad!

In 1-D math input the following four options all definitely work

#
# Three '.' characters in the range specification
# Unusual and I normally wouldn't do it but in 1-D
# it works in about 10 sec. The problem with this
# output format is that it results in a "flat"
# sequence - so it would be awkward in any subsequent
# to determine where the factors for one entry ended
# and the factors for the next entry started
#
   seq(ifactor(n!+3), n=10...50);
#
# Keeping the three '.' characters in the range
# specification: unusual and I normally wouldn't
# do it but it works in about 10 secs. Using []
# to group the output in lists. Now the factors
# for a given entry occur in a single list, which
# is probably easier for any subsequent code to
# access
#
   seq([ifactor(n!+3)], n=10...50);
#
# The above two commands, but this time with the
# 'normal' range specifications which only require
# two '.' characters in the range specifications
#
   seq(ifactor(n!+3), n=10..50);
   seq([ifactor(n!+3)], n=10..50);

If you are using 2D-input then Carl has already explained why using three '.' characters in the range specification is a BAD idea - just don't do it! The existence of three '.' characters in the range specification means that if you choose to use 2D math input, the first two commands in the above willl fail and the second two will succeed

seq[ifactor(n!+3), n=10...50);

works for me. Takes about 20secs.

seq([[ifactor(n!+3)], n=10...50);

also works in about the same time and 'groups' the output factors better

with 64-bit Firefox 49.0.1

Still nothing - how are you doing it?

Worth a try I suppose

Nope still nothing