There is a formal paradigm that applies here called "Scales (or Levels) of Measurement". The scale determines which arithmetic and statistical operations can be performed on the measurements.

The basic (Stevens's four)^[*1] scales of measurement are (this is now often taught on day 1 in "Statistics 101"):

1. Nominal (aka Categorical): Objects being measured can be categorized into a finite number of categories, but there's no numerical relationship between the categories. Examples: {female, male}, {usable, repairable, discardable}, {too sick to be worth treating---going to die anyway; should be medically treated; not sick enough to be worth treating---going to get well anyway}. Allowed operations are set membership and set cardinality.
2. Ordinal: Objects can be placed in order on some scale, but other arithmetic operations are not appropriate. Examples: hardness scales of materials, such as MOHS; Beaufort scale of hurricane strength. Allowed operations now include <, <=, and the like.
3. Interval: Objects are measured in such a way that the differences of measurements make sense and can be compared by ratio, but the ratios of individual measurements make no sense. Examples: calendar year, Fahrenheit temperation, Celsius temperature.
4. Ratio: Measurements can be compared by ratio. Examples: the vast majority of standard measurements with units, such as years or Kelvin temperature.

Although there is much debate over the details of this paradigm, the vast majority of it is that there need to be more levels^[*1]; there's wide acceptance of the ideas that there are Levels of Measurement and that some classification of them is useful. Personally, I would always include a logarithmic level for scales such as pH, dB, Richter, etc.

When you do a units conversion in Maple, it is operating on the ratio scale. Suppose I said "The temperature of the solution increased by 18 Fahrenheit degrees. What's that in Celsius?" The answer is 10 Celsius degrees. That's not the same computation as a conversion of absolute temperatures on scales that don't have equal 0 points.

[*1] See Wikipedia article "Level of measurement".

Your method is not even close to something that would work. You need to make 2 implicitplots, one of each color.

restart:
S1:= (x,y,z)-> x^2 + y^2 + z^2 - 16:  
S2:= (x,y,z)-> (x+4)^2 + y^2 + (z-2)^2 - 49:
S1S2intersect:= solve([S1,S2](x,y,z), [x,y,z], explicit):
IntPlaneX:= unapply(rhs(S1S2intersect[1][1]), (y,z));
p1:= (x,y,z)-> `if`(x <= IntPlaneX(y,z), S1(x,y,z), undefined): 
p2:= (x,y,z)-> `if`(x > IntPlaneX(y,z), S2(x,y,z), undefined):
plots:-display(
    plots:-implicitplot3d~([p1,p2], -12..5, (-7..7)$2, grid= [50$3], color=~ [red, blue]),
    style= surface
);

If the scaling (or unit of measure) of the axes is different, then lines may not appear perpendicular even though mathematically they are. This can be adjusted by including the option scaling= constrained in the plotting command.

Please put Questions in the Questions section, not the Posts section.

Expectation E (under null hypothesis of "true randomness") is that 1/2 of runs will be length 1, 1/4 length 2, 1/8 length 3, etc.

RunLengths:= proc(L::list)
local R:= Array(1..0), r:= 1, k;
    for k from 2 to nops(L) do
        if L[k-1]=L[k] then r++ else R,= r; r:= 1 fi
    od;
    Statistics:-Tally([seq((R,= r))]);
end
:
Chi2:= proc(T::list(`=`)) #return p-value (probability of randomness)
local chi2:= 0, T0:= table(sparse, T), S:= rhs(add(T)), E:= S, s:= 0, k;
    for k while E > 10 do
        chi2+= (T0[k]-(E/= 2))^2/E;
        s+= T0[k]
    od;
    #Last category is "all remaining": 
    1-Statistics:-CDF(ChiSquare(k), chi2 + (S-s-E)^2/E)
end proc
:
#Your data:
A:= [0,1,0,1,0,1,1,0,1,1,1,0,0,1,1,1,1,0,1,0,1,0,1,0,1,0,1,1,0,1,0,1,0,
0,0,1,0,1,0,1,0,1,1,1,1,1,0,0,1,0,1,0,1,0,1,0,1,1,1,0,1,0,1,0,0,0,1,1,1,
1,0,0,0,1,0,1,0,0,0,0,0,1,1,1,0,1,0,1,0,1,0,1,0,0,0,1,1,0,1,0]
:
B:=[0,1,1,0,1,1,0,1,1,1,0,1,0,0,1,1,0,0,0,0,1,1,1,1,0,0,0,0,0,0,0,0,0,0,
0,1,1,1,0,1,0,0,0,1,0,1,0,1,1,0,1,1,0,0,0,0,0,0,1,0,1,1,1,0,0,0,1,0,1,0,
0,0,0,1,0,0,0,0,1,0,1,1,1,1,1,1,0,0,0,1,1,1,1,1,1,1,0,0,1,0]
:
evalf(Chi2(RunLengths(A)));
                        0.0005236818491

evalf(Chi2(RunLengths(B)));
                          0.6136983916

Another way is like below. I prefer this slightly over using listprocedure for two closely related reasons:

1. listprocedure gives the naive user the false impression that it's possible to compute one of these functions without computing all 4. In reality, it's simply obscuring that it's always computing all 4 together.
2. Managing your own remember table gives you finer control over the recomputation.

restart: 
ddesys := {diff(S(t),t) = -beta*S(t)*Ix(t)/N, 
diff(Ex(t),t) = beta*S(t)*Ix(t)/N - sigma*Ex(t-tau__1),
diff(Ix(t), t) = sigma*Ex(t - tau__1)- gamma*Ix(t-tau__2), 
diff(R(t),t) = gamma*Ix(t-tau__2), 
diff(Dx(t),t) = delta*Ix(t), 
S(0) = 80900, Ex(0) = 1, Ix(0) = 1, R(0) = 0, Dx(0) = 0 }:
params:= {
    beta = 4, gamma = 0.0478, sigma = 0.10, delta = 0.0005,
    N=80900, tau__1 = 1.1, tau__2 = 8.7, tau__3 = 0
}:
dsn := dsolve(eval(ddesys, params), numeric):

Dsn:= proc(t) option remember; dsn(t) end proc:

F:= t-> 
    if t::numeric then eval(op(procname)(':-t'), Dsn(t)) 
    else 'procname'(args)
    fi
:

The parts that I added are in red.

Now, you can use these numeric procedures to do things like this:

plot([F[Ex], F[Ix]], 0..9);

fsolve(t-> F[Ex](t) - 10, 0..infinity);

The definitive difference between Threads and Grid (regardless of the specific command (such as Seq) used from those packages) is that Threads uses a shared-memory environment and Grid does not. A group of processes running in parallel under Threads can all "see" the same variables at higher lexical levels. Those processes can just as well overwrite those variables, which is where potential trouble arises. Two processes overwriting the same variable at the same time will usually cause erroneous results. Those commands and procedures that are guaranteed to not do that are called "threadsafe". Many commands are actually threadsafe even though they're not formally declared as such. An example is ListTools:-Classify. Most of the high-powered symbolic commands such as solve are not threadsafe.

If it's possible to do your problem with Threads, then you should: It's faster and uses much less memory. If you cannot make your code threadsafe, but it's still parallelizable, then you must use Grid.

To upload worksheets to this forum, use the green arrow on the toolbar.
 

Download betounes2.mw

In your example 2, Maple can't solve it because the input equations do not contain _xx[2]!

In your example 3: The explicit complex number is not the same thing as the RootOf because the RootOf represents all of the roots of the polynomial.

Problem 1a.
 

Download Betounes_1a.mw

Maple is not backwards compatible all the way back to version 4.3. I think that backwards compatibility only extends back to Maple 6 (not sure about the Maple 6, but I am absolutely sure that it doesn't go back to 4.3). If you read the help page ?forall, it clearly says that the first argument should be a name or list of names.  

Any situation where a procedure body (or the right side of an arrow expression) contains indirect references to variables named the same as the parameters (or the names on the left side of the arrow) will not work. This is by design, not a bug. If the procedure body can be expressed as a single expression, you can use unapply. More complicated cases require subs.

There are three dichotomies involved here, not just the two mentioned by the other respondents, and you are confusing them:

1. Standard vs. Classic GUI interface
2. 2D Input vs. 1D input (aka Maple Input)
3. Document mode vs. Worksheet mode

The red input that you've mentioned a few times is associated with 1D input. (You can change the color if you want.) You are confusing this with the other two dichotomies.

1. Standard vs. Classic: Classic is more robust and usually delivers the output faster. This has nothing to do with the speed of the computations performed by the kernel, only with the output---already computed---being formatted for your screen. There are a huge number of commands and command options (mostly related to plotting, graphics, tabular layouts) that will not work at all in Classic. Plus, you'll be restricted to a 32-bit kernel, which seems like a major setback.

2. 2D Input vs. 1D input (aka Maple Input): The 2D Input is to me nearly unreadable (due to a poorly spaced default font) and untypeable (due to constantly shifting context and the need to use arrow keys or the mouse). It is rife with bugs. Bugs are difficult to find because they may be due to character position rather than the character itself. Automatic conversion of 2D Input to 1D plaintext form results in a unreadable "brick". Despite what the other correspondents have said, there are a huge number of great syntax innovations that were introduced in Maple 2018 and 2019 that won't work at all with 2D Input. The error messages that you get for these will tend to make no sense; it's as if the 2D parser is totally unaware of that syntax. There are a moderate number of other forms of 1D syntax that work differently and erroneously when they are tried in 2D.

3. Document vs. Worksheet: I can't say much because I've never used Document.

So, I work exclusively with Standard GUI interface, 1D input (aka Maple Input), Worksheet mode.

OrbitPartition:= module()
option
	`Conceptual author: emendes`,
	`Maple code author: Carl Love <carl.j.love@gmail.com> 2020-May-13`
;
# Naming notes:
#    Names pertaining to the 1st index have "I" in the name: Is, nIs, PairsI, AllI, etc.
#    Names pertaining to the 2nd index have "J" in the name: Js, etc.
#    S always represents a tuple_size-combination of parms.

# Option cache or option remember?
#    No procedure with an 'S' argument in this module should be given option remember because
#    it takes too much memory and thus negates all the benefits of using Iterator. If such a 
#    procedure will be called again, give it option cache.
local
	#All module-level locals other than procedures
	#---------------------------------------------
	#module locals that are copied from parameters to ModuleApply:
	nIs::posint, #number of distinct 1st indices
	degree::nonnegint, #polynomial degree
	permutations_by_index::[table, table],
	CondTables::record,

	#module locals used in condition processing:
	Is::set, #set of first indices
	SubsetsI::set(set(posint)), #nonempty proper subsets of Is
	FullDegJ::posint, #critical value for condition FullDeg

 	#other module locals:
	parms::set([posint, nonnegint]), #fundamental set of pairs of indices
	nJs::posint, #number of distinct 2nd indices 
	Combos::Iterator,
	Setup::truefalse:= false, #Has problem been initialized?

	#Exclusion conditions:
	#---------------------
	#Build a table of sets of the 2nd index of S grouped by the 1st index:
	JsbyI:= proc(S)
	option cache, threadsafe;
		subs(
			_C= (op~)~(2, ListTools:-Classify(IJ-> IJ[1], S)),
			proc(i) option cache(nIs), threadsafe; `if`(assigned(_C[i]), _C[i], {}) end proc
		)
	end proc,
		
	#Make CondTables into a remember-table procedure:
     CT:= proc(K)
	option remember;
	     subs(_T= eval(CondTables[K]), proc(k) option remember; _T[k] end proc)
	end proc,
	
	Js:= proc(S) option cache; op~(2,S) end proc,

	Condition:= table([
		"AllI"=	 (S-> op~(1,S) = Is),
		"FullDeg"= (S-> max(Js(S)) >= FullDegJ),
	   	"ValidJ"=	 (S-> CT(':-ValidJ')~({Js(S)[]}) = Is),
	   	"FullDim"= (	 
	   		S-> not ormap(II-> `union`(JsbyI(S)~(II)[]) subset CT(':-FullDim')(II), SubsetsI)
   		)   		
   	]), 
   	#Conditions to check: The order is for efficiency, but that's just based on guesses.
	Conditions:= ["AllI", "FullDeg", "ValidJ", "FullDim"], 
			
	#Decide whether a permutation's orbit is allowed to be represented:
	AllowOrbit:= S-> andmap(k-> Condition[k](S), Conditions),	
	#===============================
	#End of exclusion-condition code

	#Set permutations to identity for unmentioned indices. Form Cartesian
	#product permutation of the permutations on the individual indices.
	T:= proc(T::table, ij) option remember; `if`(assigned(T[ij]), T[ij], ij) end proc,

	#This is the extension of permutation T to a permutation of tuple-size-combinations
	#of parms.
	Symm:= S-> map(IJ-> T~(permutations_by_index, IJ), S),

	#For any tuple of lists, compute its orbit under Symm. Return a
	#representative of the orbit, if allowed.
	Orbit:= proc(S)
	local r:= S, R;
		do 
			if AllowOrbit(r) then R[r]:= () else return fi 
		until assigned(R[(r:= Symm(r))]);        
		{indices(R, 'nolist')}[1] #Return lexicographic min as representative.
	end proc,

	#Process one chunk of combinations,
	#regardless of whether using sequential or multithreaded:
	DoTask:= proc(rn::[posint,posint])  #starting Iterator rank and number of iterates
	option threadsafe;
	local
		#See extended example on help page ?Iterator,RevolvingDoorCombination. 
		Combo_:= Object(Combos, 'rank'= rn[1]),
		has:= ModuleIterator(Combo_)[1], #We don't need the 'get' procedure; only need 'has'.
		C:= output(Combo_) #the Iteraror's Array output.
	;       
		(to rn[2] while has() do Orbit(parms[[seq(C+~1)]]) od)
	end proc,
		
	ModuleApply:= proc(
		tuple_size::posint, #tuple_size=1 => module setup only
		nIs::posint,
		degree::nonnegint,
		permutations_by_index::[table, table],
		CondTables::
			record(
				ValidJ::table,
				FullDim::And(table, set([set]) &under {entries})
			),
		#whether to use sequential or multithreaded (parallel) code:
		{sequential::truefalse:= false},
		#number of task "chunks", regardless of 'sequential':
		{numtasks::posint:= kernelopts('numcpus')},
		{nocompile::truefalse:= false} #Compile the Iterator?
	)
		#initializations:
		#----------------
		#Equate module locals to this procedure's parameters of the same name, as needed:
		Setup:= false;
		thismodule:-nIs:= nIs;
		thismodule:-degree:= degree;
		thismodule:-permutations_by_index:= permutations_by_index;
		thismodule:-CondTables:= CondTables;
		
		Init(2*min(numtasks, kernelopts(numcpus)));
		if tuple_size=1 then return fi; #problem setup only
		
		Combos:= Iterator:-RevolvingDoorCombination(
			nops(parms), tuple_size, ':-compile'= not nocompile
		);
		
		#Iterations:
		for local iter to 2 do #kludge to workaround error on 1st run after restart
			try
				return
					#Iterate over all "chunks" of combinations. (Chunks are
					#task-sized subsets of combinations.)
					`if`(sequential, map, Threads:-Map['tasksize'= 1])(
						DoTask,			
						{Iterator:-SplitRanks(Number(Combos), ':-numtasks'= numtasks)[]}
					)
			catch:
				if iter=1 then
					printf(
						"Error trapped, "
						"so efficiency measures are invalid for this run.\n"
					)
				else
					error
				fi
			end try
		od
	end proc,

	Init:= proc(cachesize::posint)
	local P, nterms:= (nIs,deg)-> binomial(nIs+deg, deg);
		Is:= {$1..nIs}; #1st indices
		FullDegJ:= nterms(nIs, degree-1);
		local Js:= {$0...(nJs:= nterms(nIs, degree))-1}, i, j;
		parms:= {seq(seq([i,j], i= Is), j= Js)};
		SubsetsI:= {seq(combinat:-choose(Is, i)[], i= 1..nIs-1)};
		forget~([for P in thismodule do P od])[];
		for P in [thismodule:-Js, JsbyI] do
			try 
				Cache:-Resize(cachesize, P)
			catch "procedure should have a cache":
				Cache(cachesize, procedure= P)
			end try
		od;		
		
		#Validate permutations:
		local ind, ent;
		for local ij,PT in permutations_by_index do
			if
				nops((ind:= {indices(PT, 'nolist')})) <>
				nops((ent:= {entries(PT, 'nolist')}))
				or not ind union ent subset op~(ij, parms)
			then
				printf("Index %d = %a.", ij, [Is,Js][ij]);
				error "%-1 table is not a permutation of its index", ij
			fi
		od;
		
		#Validate CondTables:
		if {indices(CondTables:-ValidJ, 'nolist')} <> Js then
			error "indices of ValidJ must be %1", Js
		fi;
		if not {entries(CondTables:-ValidJ, 'nolist')} subset Is then
			error "entries of ValidJ must be *sequences* from %1", Is
		fi;
		if {indices(CondTables:-FullDim, 'nolist')} <> SubsetsI then
			error "indices of FullDim must be %1", PairsI
		fi;
		if not `union`(entries(CondTables:-FullDim, 'nolist')) subset Js then
			error "entries of FullDim must be subsets of %1", Js
		fi;

		Setup:= true;
		return
	end proc,

	ModuleLoad:= proc() Setup:= false; return end proc
;
export
	CondCheck:= proc(S::set([posint, nonnegint]))
		if not Setup then error "first, call OrbitPartition to setup problem" fi;
		if not S::set([integer[1..nIs], integer[0..nJs-1]]) then
			error "index out range; nIs=%1, nJs=%2", nIs, nJs 
		fi;
		map(evalb@apply, Condition, S)
	end proc
;
	ModuleLoad()	
end module:

Download Symmetries.mw

A completely different solution is to extract the "y" values from a plot of the function. You don't need to look at the plot.

seq(op([1,1], plot(f, 0..1))[..,2]);

For finer control, this can be used with plot options numpoints, adaptive, and\or sample (see ?plot,options).

Plots can be combined with plots:-display. In the plot below, I used y1 > 0 for the first plot and y2 < 0 for the second so that they don't overlap.

restart:
initialset1 := {seq(seq([x1(0) = a, y1(0) = b], a = -2 .. 2), b = -2 .. 2)}:
DE1 := [diff(x1(t), t) = y1(t), diff(y1(t), t) = 1-x1(t)]:
p1:= DEtools:-phaseportrait(
    DE1, [x1, y1], t = -5 .. 5, initialset1, x1 = -2 .. 2, y1 = 0 .. 3, 
    stepsize = .1, color = green, numpoints = 600, thickness = 2, 
    linecolor = black
):
initialset2 := {seq(seq([x1(0) = a, y1(0) = b], a = -2 .. 2), b = -2 .. 2)}:
DE2 := [diff(x1(t), t) = -y1(t), diff(y1(t), t) = 1+x1(t)]:
p2:= DEtools:-phaseportrait(
    DE2, [x1, y1], t = -5 .. 5, initialset2, x1 = -2 .. 2, y1 = -3 .. 0, 
    stepsize = .1, color = green, numpoints = 600, thickness = 2, 
    linecolor = black
):
plots:-display(p1,p2);

There: That's substantial progress towards your goal. Now you just need to adjust the initial conditions so that the trajectories connect across the x-axis.